Higham - Accuracy and Stability of Numerical Algorithms (523152), страница 102
Текст из файла (страница 102)
ThenButThereforeThere is still the question ofwhether Xn is optimally scaled. We can perturb the zero eigenvalues of A to distinct,sufficiently small numbers so that ||Ak||2 is essentially unchanged and so that theonly freedom in Xn is a diagonal scaling. Since Xn has columns of unit 2-norm,rninF=diag(fi) k2(XnF) > n–1/2 k2(Xn) (Theorem 7.5), so even for the optimallyscaled Xn the bound can be arbitrarily poor.17.2. A simple example isfor which λ(A) = {–1,1,3} and λ(|A|) = {–l,2±(4α+l)1/2}, so p(|A|)/p(A)S OLUTIONSTOP ROBLEMS56918.1. For the Householder matrix I – (2/vTv)vvT there is an eigenvalue –1 witheigenvector v, and there are n – 1 eigenvalues 1, with eigenvectors any basis forspan(v).A Givens matrix G(i,j,θ) has eigenvalues e± iθ, together with n – 2eigenvalues 1.18.2. Straightforward manipulation shows that a bound holds of the form cnu +0(u2), where c is a constant of order 10.18.3.
We must have x*Px = x*y, so x*y must be real. Also, we need x*x = y*yand x y. Then, with v = x – y, v*v = 2v*x, so Px = y.LAPACK uses modified Householder matrices P = I – tvv* that are unitarybut not Hermitian (T is not real). The benefit of this approach is that such a Pcan always be chosen so that Px = αe1, with α real (for this equation to hold fora genuine Householder matrix P it would be necessary that x1 be real). For moredetails see Lehoucq [698, 1994].18.4. False.
det(P) = – 1 for any Householder matrix and det(G) = 1 for anyGivens matrix, so det(Q) = 1. Moreover, whereas P is generally a full, symmetricmatrix, Q has some zero entries and is generally nonsymmetric.18.5. The inequalities follow from the fact that, in the notation of (18.1), ||xk||2 =|rkk| and ||ck(:,j)|| =the latter equality being a consequence of theinvariance of the 2-norm under orthogonal transformations. The last part followsbecause QR factorization with column pivoting on A is essentially equivalent toCholesky factorization with complete pivoting applied to ATA.18.6.
If y = |W|x and W comprises r disjoint rotations then, in a rather loosenot at ion,18.7. Straightforward. This problem shows that the CGS and MGS methodscorrespond to two different ways of representing the orthogonal projection ontospan{q1, . . . ,qj}.18.8. Assume, without loss of generality, that ||a1||2 < ||a2||2. If E is any matrixsuch that A + E is rank deficient thenWe take E = [el, 0], where el is chosen so that= 0 and al + el = αa2, forsome α.
From Pythagoras’s theorem we have that ||el||2 = tan θ ||al||2, and soTogether with the trivial boundthe result.this yieldsS OLUTIONS570TOP ROBLEMS18.9. We find thatfor i j, but forIt is easy to see= 1 impliesHence forand thatshowing, as expected, that the loss of orthogonality for MGSis bounded in terms of k2(A)uForthat18.10. P is the product Pn . . . P1, where Pk is defined in (18.20). Since= 0, we haveHenceThis is of the required form, with Q = [q1,.. ., qn] the matrix from the MGS methodapplied to A.18.11. With Q defined as in the hint,∆ A=====QR–A=(Q–P21)R+∆A 2T(VW – P21)R + ∆A2V(I – S)WTR + ∆AzV(I + S)-1 C2 WTR + ∆A2V(I + S)-l WT · WCUT · UCWT · R + ΑAzand the result follows.18.12. To produce the Householder method’s P we have to explicitly form theproduct of the individual Householder transformations.
As long as this is done ina stable way the computed P is guaranteed to be nearly orthogonal. MGS’s Q isformed in an algebraically different way, and the rounding errors in its formation aredifferent from those in the formation of P; in other words, Q is not a submatrix ofP. Consequently, there is no reason to expect Q to be nearly orthonormal. Furtherinsight can be gained by a detailed look at the structure of the computed P; seeBjörck and Paige [119, 1992, §4].S OLUTIONSTOP ROBLEMS57118.13. It is straightforward to show that ATA– I = (A – U)T(A + U).
Takingnorms gives the lower bound. Since A + U = U(H + I) we have, from the previousrelation,(A - U)TU = (ATA - I)(H + I)-1.HenceIn fact, the result holds for any unitarily invariant norm (but the ||A||2 + 1 in thedenominator must be retained).19.1. One approach is to let x be a solution and y an arbitrary vector, and considerf (α) = ||A(x + αy) – b|| By expanding this expression it can be seen that if thenormal equations are not satisfied then α and y can be chosen so that f(α) < f(0).Alternatively, note that for f(z) = (b – Ax)T(b – Ax) = xTATAx – 2bTAx + bTbwe have f(x) = 2(ATAx – ATb) and 2f(x) = 2ATA, so any solution of thenormal equations is a local minimum of f, and hence a global minimum since f isquadratic.
The normal equations can be written as (b – Ax)TA = 0, which showsthat the residual b – Ax is orthogonal to range(A).19.2. A = Q [~]. The computed upper triangular QR factor satisfies A +∆A =with |∆A| < mnγ cmG1|Al and ||G1||F = 1, by Theorem 18.4. BYLemma 18.3 the computed transformed right-hand side satisfies= QT(b + ∆b),with |∆b| < mnγ cmG1 |b|.By Theorem 8.5, the computed solutionto the triangular systemsatisfiesSoexactly solves the LS problemwhich is equivalent to the LS problem, on premultiplying by Q,We havewhere G > G1 and ||G||F = 1. The normwise bounds are proved similarly.19.3.
A straightforward verification.S OLUTIONS572TOP ROBLEMS19.4. Notice thatwhich is minimized if ||Axi — ei ||2 is minimized for i = 1: m. Thus we have m+independent LS problems. The solution is xi = A ei , i = 1: m, that is, X =++A Im = A . This is the unique solution only if A has full rank, but it is always theunique minimum Frobenius norm solution.19.5. By Theorem 18.12 there is an orthonormal matrix [W1, wn+l]such thatwithThe computed solutionsatisfiesThereforeexactly solves the LS problemwhich, on premultiplying by [W1, wn+1], is equivalent to the LS problemWe havewhereand= 1. Normwise bounds on19.6. Subtractingfromandfollow similarly.we haveby (19.12).Taking norms gives the result.toS OLUTIONSTOP ROBLEMS57319.7.
By constructing a block LU factorization it is straightforward to show thatdet(C(α) – λI) = (α – λ)m-n det(ATA + λ(α – λ)).Hence the eigenvalues of C(α) are λ = α (m – n times) together with the solutionsof λ(α – λ) =namely,Nowso the minimum condition number must occur whengivesfor whichThisThe lower bound for the maximum is achieved for19.8. Let y = 0 and ||(A + ∆A)y – b||2 = min. If b = 0 then we can takeAA= O. Otherwise, the normal equations tell us that (A+ ∆A)Tb = 0, so ||∆ A||2 >||ATb||2/||b||2. This lower bound is achieved, and the normal equations satisfied, for∆AT = –(ATb)bT/bTb. Hence19.9. For the case λ∗ < 0 we have20.1.
Setting the gradient of (20.13) to zero gives ATAx – ATb + c = 0, whichcan be written as y = b – Ax, ATy = c, which is (20.12). The Lagrangian for(20.14) is L(y, x) = ½(y – b)T(y – b) + xT(ATy – c).yL(y, x) = y – b + Ax, andTxL(y, x) = A y – c. Setting the gradients to zero gives the system (20.12).21.1. The modified version has the same flop count as the original version.21.4. The summation givesS OLUTIONS574TOP ROBLEMSwhich implies the desired equality.
It follows that all columns of V–1 except thefirst must have both positive and negative entries. In particular, V-1 > 0 is notpossible. The elements of V-1 sum to 1, independently of the points αi (see alsoProblem 13.7).21.5. We have U(i,:)T(αo,. . . ,αn) = W(i,:)V(αo,. . . ,αn)- But T = LV, whereL has the form illustrated for n = 5 byT-TT-Tand Le = e, so U(i,:)L = W(i,:), or U(i,:) = L W(i,:) .
But L > 0-T-1-lnby the given x formula, so ||L ||1 = ||L || = ||L e|| = ||e|| = 1, hence||U(i,:)||1 < ||W(i,:)||1.As an aside, we can evaluate ||L|| asafter a little trigonometric algebra.21.7. Denote the matrix by C. For the zeros of Tn+1 we haveIt follows thatExtrema of Tn: we havegiving the result.CDCT =Hence B = CD½ satisfies BBT =k2 (C) < k2(D½)k2 (B) = 2.and so k 2 (B) = K 2 (D) ½ =Then21.8. The increasing ordering is never produced, since the algorithm must chooseα1 to maximize |α l – α0|.21.10.
The dual condition number isSee Higham [533, 1987] for proofs.S OLUTIONSTOP ROBLEMS57522.1. Consider the case where m < min(n,p), and suppose n = jm and p = km forsome integers j and k. Then the multiplication of the two matrices can be split intom x m blocks:which involves a total of jk multiplications of m x m matrices, each involving O(mα )operations.
Thus the total number of operations is O(jkmα ), or O(mα -2np), asrequired, and we can show similar results for the cases when n and p are the smallestdimensions.22.2. ½n3 + n2 multiplications and ³/2n3 + 2n(n – 1) additions.22.3. For large n = 2k, Sn (8)/Sn (n) 1.96 x (7/8) k and S n (l)/S n (8) 1.79.The ratio Sn(8)/Sn(n) measures the best possible reduction in the amount of arithmetic by using Strassen’s method in place of conventional multiplication. The ratioS n (1)/S n (8) measures how much more arithmetic is used by recurring down to thescalar level instead of stopping once the optimal dimension n0 is reached.
Of course,the efficiency of a practical code for Strassen’s method also depends on the variousnon-floating-point operations.22.5. Apart from the differences in stability, the key difference is that Winograd’sformula relies on commutativity and so cannot be generalized to matrices.22.7. Some brief comments are given by Douglas, Heroux, Slishman, and Smith[317, 1994].22.9. The inverse isHence we can form AB by inverting a matrix of thrice the dimension. This resultis not of practical value, but it is useful for computational complexity analysis.24.2.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
