William W. Hager. Condition estimates. SIAM Journal on Sci. Statist. Comput., 5(2):311-- 316, 1984.  Home/Search   Document Not in Database   Summary   Related Articles   Check

.... M N) by solving triangular matrix equations: kxk 2 kck 2 kXkF kCkF kZ SYCT k 2 = min (ZSYCT ) Sep The right hand side C is chosen such that the lower bound gets as large as possible. This leads to a Frobenius norm based estimate. For computation of 1 norm based estimates see [17, 19, 31]. The Sep functions associated with the Sylvester type matrix equations are: Sep[SYCT] inf kXkF =1 kAX XBkF = min (ZSYCT ) Sep[LYCT] inf kXkF =1 kAX X( A )k F = min (ZLYCT ) Sep[GCSY] inf k(X;Y )kF =1 k(AX Y B; DX Y E)kF = min (ZGCSY ) The same techniques as presented above ....

W.W. Hager. Condition Estimates, SIAM J. Sci. Stat. Comp., 5:311-316, 1984.

....for the order of magnitude of k(A) 30] After such a proposal, several researchers have attempted to enhance the previous result [31, 32] until, in 1984, W. Hager suggested a new technique for estimating the l 1 condition number with an error under 5 = oo and a probability higher than 0. 97 [33]. Finally, it is worth mentioning that, because of the very nature of the approximation method, all such algorithms calculate an under estimation of the actual condition number. The algorithm can easily be extended to the general case of nonlinear equations. 13 A solution x such that E( x) ....

W. W. Hager, "Condition estimates," SIAM J. Sci. Stat. Comput, vol. 5, no. 2, pp. 311--316, 1984.

....number of variables, one ends up at some stage with a sparse linear system which has at least interval coe cients on the right hand side. Or, in noninterval terms, one gets expressions involving the absolute value or norm of the inverse. 3 There are techniques based on an observation of Hager [6] (see also Neumaier [19, Algorithm 2.5.5] that frequently produce good approximations or even exact values for the norm of the inverse. But as with sensitivity analysis or Monte Carlo methods, one never knows when this is the case, and occationally Hager type estimators fail by orders of ....

W.W. Hager, Condition estimates, SIAM J. Sci. Statist. Comput. 5 (1984), 311-316.

....null vector (namely, A x) and it is this rather than the kA k estimate that is required in some applications. Indeed a precursor of the LINPACK algorithm was developed for such an application [71] The now widespread use of condition estimation relies on two developments. First, Hager [76] devised a method for estimating kBk 1 that computes a (usually small) number of matrix vector products involving B and B . An immediate advantage over the LINPACK estimator is that Hager s method can be programmed as a black box that requires no knowledge of the details of how the products Bx ....

William W. Hager. Condition estimates. SIAM J. Sci. Stat. Comput., 5 (2):311--316, 1984.

....(function rcond) During the 1980s, attention was drawn to various componentwise condition numbers and it was recognized that most condition estimation problems can be reduced to the estimation of kAk when matrix vector products Ax and A x can be cheaply computed [2] 16, Sec. 14.1] Hager [12] derived an algorithm for the 1 norm that is a special case of the more general p norm power method proposed by Boyd [3] and later investigated by Tao [18] Hager s algorithm was modified by Higham [14] and incorporated in LAPACK (routine xLACON) 1] and Matlab (function condest) The LINPACK and ....

....we note that although our work is specific to the 1 norm, the 1 norm can be estimated by applying our algorithm to A , since kAk1 = kA k 1 . 2. Block 1 Norm Power Method. The 1 norm power method is a special case of Boyd s p norm power method [3] and was derived independently by Hager [12]. For a real matrix A we denote by sign(A) the matrix with (i; j) element 1 or Gamma1 according as a ij 0 or a ij 0. The jth column of the identity matrix is denoted by e j . Algorithm 2.1 (1 norm power method) Given A 2 R this algorithm computes fl and x such that fl kAk 1 and kAxk 1 = ....

[Article contains additional citation context not shown here]

William W. Hager. Condition estimates. SIAM J. Sci. Stat. Comput., 5(2):311--316, 1984.

.... jbj) k1 : 5. 1) The numerator in the bound is of the form k jA jd k1 , and as in [2] we have jd k1 = k jA jDe k1 = k jA Dje k1 = k jA Dj k1 = k A D k1 ; where D = diag(d) and e = 1; 1; 1) Hence k jA jd k1 can be estimated using the norm estimator of [11, 18, 19], which estimates kBk 1 at the cost of forming a few matrix vector products involving B and B . With B = A D) we need to solve a few linear systems involving A and A . The bound (5.1) is the one returned by the linear equation solvers in the Fortran linear algebra library LAPACK [1] ....

....decomposition, using the level 3 BLAS routine xGEMM to transform the right hand side C, calling xTRSYL to solve the (quasi ) triangular Sylvester equation, and using xGEMM to transform back to the solution X . The error bound (5. 2) can be estimated using xLACON (which implements the estimator of [11, 18, 19]) in conjunction with the above routines. We have written a Fortran 77 code dggsvx that follows the above outline. It is in the style of an LAPACK driver and follows the LAPACK naming conventions. Acknowledgements. I thank Zhaojun Bai for bringing the question of backward error for the Sylvester ....

William W. Hager, Condition estimates, SIAM J. Sci. Stat. Comput., 5 (1984), pp. 311--316.

....an important tool for assessing the quality of the computed solutions. Typically, these condition numbers are as expensive to compute as the solution itself [8] The LAPACK [1] and ScaLAPACK [2] condition numbers and error bounds are based on estimated condition numbers, using the method of Hager [5], which was subsequently improved by Higham [6] Hager s method estimates kBk 1 given only the ability to compute matrix vector products Bx and B T y. If we take B = A 1 and compute the required products by solving linear systems with A, we obtain an estimate of the 1 norm condition number ....

....numbers can be arbitrarily poor estimates exist [6] 7] Moreover, when the accuracy of the estimates becomes important for certain applications [9] the method does not provide an obvious way to improve the estimate. Higham and Tisseur [9] present a block generalization of the estimator of [5], 6] that iterates with an n t matrix, where t 1 is a parameter, enabling the exploitation of matrix matrix operations (level 3 BLAS) and thus promising greater eciency and parallelism. The block algorithm also o ers the potential of better estimates and a faster convergence rate, through ....

W. W. Hager. Conditions estimates. SIAM J. Sci. Stat. Comput., 5:311-316, 1984.

....an important tool for assessing the quality of the computed solutions. Typically, these condition numbers are as expensive to compute as the solution itself [6] The LAPACK [1] and ScaLAPACK [2] condition numbers and error bounds are based on estimated condition numbers, using the method of Hager [3], which was subsequently improved by Higham [4] Hager s method estimates kBk 1 given only the ability to compute matrix vector products Bx and B T y. If we take B = A 1 and compute the required products by solving linear systems with A, we obtain an estimate of the 1 norm condition number 1 ....

....numbers can be arbitrarily poor estimates exist [4] 5] Moreover, when the accuracy of the estimates becomes important for certain applications [7] the method does not provide an obvious way to improve the estimate. Higham and Tisseur [7] present a block generalization of the estimator of [3], 4] that iterates with an n t matrix, where t 1 is a parameter, enabling the exploitation of matrix matrix operations (level 3 BLAS) and thus promising greater eciency and parallelism. The block algorithm also o ers the potential of better estimates and a faster convergence rate, through ....

W. W. Hager. Conditions estimates. SIAM J. Sci. Stat. Comput., 5:311-316, 1984.

.... T jjA T xj and G = diag(g i ) we have k jI Gamma A Aj Delta jA T j Delta jA T xj k1 = k jI Gamma A Ajg k1 = k jI Gamma A AjGe k1 = k jI Gamma A AjG k1 = k j(I Gamma A A)Gj k1 = k (I Gamma A A)G k1 : The latter norm can be estimated by the method of [8] and [9, 11] which estimates kBk 1 given a means for forming matrix vector products Bx and B T y. Forming these products for B T = I Gamma A A)G involves multiplying by G and Q, or their transposes, and solving triangular systems with R and R T . 10 Table 4.1 A = randsvd( 10; 16] ....

W.W. Hager, Condition estimates, SIAM J. Sci. Statist. Comput., 5 (1984), pp. 311--316.

.... i ) It is easy to show that 1 2 1 (A; x) cond1 (A; x) 1 (A; x) The quantity 1 (A; x) can be estimated without explicitly forming the matrices A Gamma1 XBD 1 2 IR n Thetat and A Gamma1 D 2 2 IR n Thetan (assuming a factorization of A is available) by using the method of Hager [10] and Higham [12, 14] this method estimates kCk1 at the cost of forming a few matrix vector products Cx and C T y. We mention also two interesting nonlinear structures, those of Vandermonde matrices V = ff i Gamma1 j ) and Cauchy matrices H = Gamma (ff i fi j ) Gamma1 Delta . In ....

William W. Hager, Condition estimates, SIAM J. Sci. Stat. Comput., 5 (1984), pp. 311--316.

....y = Ax z = A T dual p (y) if kzk q z T x fl = kyk p quit end x = dual q (z) end Algorithm PM requires about 4rmn flops if there are r iterations for convergence. The convergence test can be written in several different ways, as we explain below; the form chosen here is the one used in [10, 13, 14]. The power method was first derived and analysed by Boyd [3] and it was later investigated by Tao [26] Tao applies the method to an arbitrary mixed subordinate norm kAk ff;fi = max x6=0 kAxk ff kxk fi ; 2.1) while Boyd takes the ff and fi norms to be p norms (possibly different) Algorithm ....

....with fl = kAk p = kxk p kyk q , whatever x 0 . 2) Boyd [3] shows that if A has nonnegative elements, A T A is irreducible, 1 p 1, and x 0 has positive elements, then the x k converge and fl k kAk p . In the case p = 1, Algorithm PM is a 1 norm estimation algorithm devised by Hager [10] (independently of [3] and [26] and subsequently analysed and modified by the present author [13, 14] The algorithm given in [13] is the basis of all condition number estimation in LAPACK [2] Algorithm PM has two remarkable properties when p = 1: it almost always converges within four ....

W.W. Hager, Condition estimates, SIAM J. Sci. Stat. Comput., 5 (1984), pp. 311-- 316.

.... of these counter examples the LINPACK estimator has been widely used and is regarded as being almost certain to produce an estimate correct to within a factor ten in practice [27] Another 1 norm condition estimation algorithm was developed by Higham [29, 30] building on an algorithm of Hager [25]. This estimator is in the NAG library and is being used throughout LAPACK [2] The general algorithm estimates kBk 1 given a means for forming matrix vector products Bx and B T y. By taking B = A Gamma1 and using an LU factorization of A we obtain an estimator with the same functionality as ....

William W. Hager, Condition estimates, SIAM J. Sci. Stat. Comput., 5 (1984), pp. 311--316.

....j kxk Gamma1 1 I and we recover the expression for cond 1 (A; x) in (1.2) Finally, we note that when p = 1 or 1, the bounds in (3. 12) can be estimated in O(n 2 r) operations without forming A Gamma1 if a QR or LU factorization of A is available; this can be done using the method of Hager [15] and Higham [17, 18] which estimates kBk 1 or kBk1 by evaluating several matrix vector products involving B and B T . The use of this method to estimate a componentwise condition number was first suggested in [2] in connection with the condition number (1.2) and the latter condition number is ....

W.W. Hager, Condition estimates, SIAM J. Sci. Statist. Comput., 5 (1984), pp. 311--316.

....Demmel and Duff [2] In particular, they show how to estimate the condition number (2.10) cheaply, given a factorization of A. Their approach avoids explicit computation of A Gamma1 by manipulating the condition number into a form that can be estimated using a matrix norm estimator of Hager [36] and Higham [40] This approach can be adapted to estimate virtually any form of componentwise condition number for a linear system. Finally, we mention how the analysis of this section is reflected in LAPACK. The LAPACK expert driver routine xGESVX 1 solves a linear system Ax = b by 1 The ....

William W. Hager, Condition estimates, SIAM J. Sci. Stat. Comput. 5 (1984), no. 2, 311--316.

....developed see [21] for a survey. The most well known estimator is the one used in LINPACK, which provides a lower bound for 1 (A) The method underlying this estimator does not generalize to the estimation of E;f (A; b) A more versatile estimator with this capability is one developed by Hager [20] and Higham [23] This estimator treats the general problem of estimating kBk 1 , where B is not known explicitly. The estimator assumes that B is described by a black box that can evaluate Bx or B T x given x. Typically, 4 or 5 such matrix vector products are required to produce a lower bound ....

W.W. Hager, Condition estimates, SIAM J. Sci. Stat. Comput., 5 (1984), pp. 311--316.

.... GammaB is very small, then the change in X can be very large. Note that sep(A; GammaB) oe min (T) 1=kT Gamma1 k 2 where oe min (T) is the minimum singular value of T = I n Omega A B 0 Omega I m . The quantity sep(A; GammaB) is efficiently estimated by the algorithm described in [6], 8] which is implemented in LAPACK. For this purpose a few linear systems with coefficient matrices T and T 0 are solved, which is equivalent to solve a few Sylvester equations AX XB = C and A 0 X XB 0 = C. Using the fact that [11] sep(A; GammaB) sep(U 0 AU; GammaV 0 BV) it is ....

W. W. Hager. Condition estimates. SIAM J. Sci. Stat. Comput., 5:311-- 316, 1984.

....iterative and require estimates of the condition number of some triangular submatrices at every iteration step of initial factorization, refinement and updating. In the URV and the ULV factorizations small singular values and associated null vectors are estimated by means of conditions estimators [3, 4, 5, 15, 24, 30]. A survey of condition estimators is given in [16] In the practical ULV (resp. URV ) factorization, however, each left (resp. right) null vector is recomputed from its corresponding right (resp. left) null vector via triangular solves. Triangular solves are required for the initial ....

W. W. Hager, Condition estimates, SIAM J. Sci. Statist. Comput. 29, (1987), 311-316.

....(function rcond) During the 1980s, attention was drawn to various componentwise condition numbers and it was recognized that most condition estimation problems can be reduced to the estimation of kAk when matrix vector products Ax and A T x can be cheaply computed [2] 16, Sec. 14.1] Hager [12] derived an algorithm for the 1 norm that is a special case of the more general p norm power method proposed by Boyd [3] and later investigated by Tao [18] Hager s algorithm was modified by Higham [14] and incorporated in LAPACK (routine xLACON) 1] and Matlab (function condest) The LINPACK and ....

....we note that although our work is specific to the 1 norm, the 1 norm can be estimated by applying our algorithm to A , since kAk1 = kA k 1 . 2. Block 1 Norm Power Method. The 1 norm power method is a special case of Boyd s p norm power method [3] and was derived independently by Hager [12]. For a real matrix A we denote by sign(A) the matrix with (i; j) element 1 or Gamma1 according as a ij 0 or a ij 0. The jth column of the identity matrix is denoted by e j . Algorithm 2.1 (1 norm power method) Given A 2 R n Thetan this algorithm computes fl and x such that fl kAk 1 and ....

[Article contains additional citation context not shown here]

William W. Hager. Condition estimates. SIAM J. Sci. Stat. Comput., 5(2):311--316, 1984.

....of the form k jA Gamma1 jd k1 , and as in [2] we have k jA Gamma1 jd k1 = k jA Gamma1 jDe k1 = k jA Gamma1 Dje k1 = k jA Gamma1 Dj k1 = k A Gamma1 D k1 ; where D = diag(d) and e = 1; 1; 1) T . Hence k jA Gamma1 jd k1 can be estimated using the norm estimator of [11, 18, 19], which estimates kBk 1 at the cost of forming a few matrix vector products involving B and B T . With B = A Gamma1 D) T we need to solve a few linear systems involving A and A T . The bound (5.1) is the one returned by the linear equation solvers in the Fortran linear algebra library ....

....decomposition, using the level 3 BLAS routine xGEMM to transform the right hand side C, calling xTRSYL to solve the (quasi ) triangular Sylvester equation, and using xGEMM to transform back to the solution X . The error bound (5. 2) can be estimated using xLACON (which implements the estimator of [11, 18, 19]) in conjunction with the above routines. We have written a Fortran 77 code dggsvx that follows the above outline. It is in the style of an LAPACK driver and follows the LAPACK naming conventions. Acknowledgements. I thank Zhaojun Bai for bringing the question of backward error for the Sylvester ....

William W. Hager, Condition estimates, SIAM J. Sci. Stat. Comput., 5 (1984), pp. 311--316.

.... k xk 2 = hm 1m je T m yj kAk F k xk 2 ; and we do not have to compute the residual kA x Gamma xk 2 explicitly as long as j AE u: The cost of computing kAk F is significant for large matrices and an interesting way to reduce this is to compute an estimate of this norm as developed in Hager (1984) and Higham (1990) A rough estimate of kAk can often be obtained from physical considerations when A is the discretization matrix of some partial differential equation. Acknowledgments We thank Prof. Chatelin (University of Paris IX and CERFACS) Prof. I. Duff (CERFACS and RAL) Dr. J. Scott ....

W. W. Hager, (1984), Condition estimates, SIAM J. Sci. Stat. Comput., 5, 311-- 316.

....is of the form k jA Gamma1 jd k1 , and as in [2] we have k jA Gamma1 jd k1 = k jA Gamma1 jDe k1 = k jA Gamma1 Dje k1 = k jA Gamma1 Dj k1 = k A Gamma1 D k1 ; where D = diag(d) and e = 1; 1; 1) T . Hence k jA Gamma1 jd k1 can be estimated using the norm estimator of [11, 18, 19], which estimates kBk 1 at the cost of forming a few matrix vector products involving B and B T . With B = A Gamma1 D) T we need to solve a few linear systems involving A and A T . The bound (5.1) is the one returned by the linear equation solvers in the Fortran linear algebra library ....

....decomposition, using the level 3 BLAS routine xGEMM to transform the right hand side C, calling xTRSYL to solve the (quasi ) triangular Sylvester equation, and using xGEMM to transform back to the solution X. The error bound (5. 2) can be estimated using xLACON (which implements the estimator of [11, 18, 19]) in conjunction with the above routines. We have written a Fortran 77 code dggsvx that follows the above outline. It is in the style of an LAPACK driver and follows the LAPACK naming conventions. Its leading comment lines together with an example program are listed in the appendix. ....

W.W. Hager, Condition estimates, SIAM J. Sci. Statist. Comput., 5 (1984), pp. 311--316.

....factorization of A. We discuss situations in which scaling is necessary to prevent overflow and give an example of how our routines are used in the LAPACK condition estimators. 1 Introduction The condition estimators in the LAPACK software package are based on the estimate described by Hager [3] and Higham [4] The norm of the matrix A is computed by conventional means, an estimate is obtained for the 1 norm (or infinitynorm) of its inverse, and the reciprocal condition estimate is computed as RCOND = 1= kAkkA Gamma1 k) When used to estimate kA Gamma1 k, Higham s version of ....

....routines. 2 Hager s algorithm for estimating the norm of a matrix In this section, we briefly describe Hager s method for estimating the 1 norm of a matrix. The main application is in estimating the norm of B = A Gamma1 when B has not been explicitly computed. For further details, refer to [3] or [4] The 1 norm of an n by n matrix B is defined in terms of the vector 1 norm kxk 1 = P n i=1 jx i j as kBk 1 = max x6=0 kBxk 1 kxk 1 : A well known property is that kBk 1 is equal to the maximum column sum: kBk 1 = max 1jn n X i=1 jb i;j j; 2:1) hence the maximum in (2.1) is ....

W. W. Hager. Condition estimates. SIAM J. Sci. Stat. Comput., 5:311--316, 1984.

.... T xj and G = diag(g i ) we have k jI Gamma A Aj Delta jA T j Delta jA T xj k 1 = k jI Gamma A Ajg k 1 = k jI Gamma A AjGe k 1 = k jI Gamma A AjG k 1 = k j(I Gamma A A)Gj k 1 = k (I Gamma A A)G k 1 : The latter norm can be estimated by the method of [8] and [9, 11] which estimates kBk 1 given a means for forming matrix vector products Bx and B T y. Forming these products for B T = I Gamma A A)G involves multiplying by G and Q or their transposes, and solving triangular systems with R and R T . Table 4.1: A = randsvd( 10; 16] ....

W.W. Hager, Condition estimates, SIAM J. Sci. Statist. Comput., 5 (1984), pp. 311--316.

No context found.

W. W. Hager, Condition estimates, SIAM J. Sci. Statist. Comput., 5 (1984), pp. 311--316.

No context found.

W. W. Hager, Condition estimates, SIAM J. Sci. Comput., 5 (1984), pp. 311-316.

....For the airline scheduling matrix of section 7, L has up to 1.49 million nonzeros and it is impractical to compute the product L D L T after each update. To obtain a quick and accurate estimate for #E# 1 , where E = AA T L D L T , we applied the strategy presented in [23] (see [24, p. 139] for a symbolic statement of the algorithm) to estimate the 1 norm of a matrix. That is, we used a gradient ascent approach to compute a local maximum for the following problem: max #Ex# 1 : #x# 1 = 1 . Since L is used multiple times in the following algorithm, we copied our ....

....8000 10000 12000 14000 0 0.5 1 1.5 2 2.5 3 3.5 4 x 10 12 step Fig. 7.1. Estimated 1 norm of error in the LDL T factorization. To improve the accuracy of the 1 norm estimate, we used Algorithm 7 three times. In the second and third trials, a di#erent starting vector x was used as described in [23]. Observe that Algorithm 7 only makes use of the product between the matrix E and a vector. This feature is important in the context of sparse matrices since E contains the term L D L T . It is impractical to compute the product L D L T , but it is practical to multiply L D L T by a ....

W. W. Hager, Condition estimates, SIAM J. Sci. Stat. Comput., 5 (1984), pp. 311--316.

....factorization. For the airline scheduling matrix of x8, L has up to 1.49 million nonzeros and it is impractical to compute the product L L T after each update. To obtain a quick and accurate estimate for kEk 1 , where E = AA T Gamma L L T , we applied the strategy presented in [16] (see [17, p. 139] for a symbolic statement of the algorithm) to estimate the 1 norm of the inverse of a matrix. That is, we used a gradient accent approach to compute a local maximum for the following problem: maxfkExk 1 : kxk 1 = 1g: Since L is used multiple times in the following algorithm, ....

....E T y j = arg max fjz i j : i = 1 to mg if jz j j z T x return x i = 0 for i = 1 to n x j = 1 end while end Algorithm 8 To improve the accuracy of the 1 norm estimate, we used Algorithm 8 three times. In the second and third trials, a different starting vector x was used as described in [16]. Observe that Algorithm 8 only makes use of the product between the matrix E and a vector. This feature is important in the context of sparse matrices since E contains the term L L T and it is impractical to compute the product L L T , but it is practical to multiply L L T by a ....

W. W. Hager, Condition estimates, SIAM J. Sci. Comput., 5 (1984), pp. 311-- 316.

No context found.

William W. Hager. Condition estimates. SIAM Journal on Sci. Statist. Comput., 5(2):311-- 316, 1984.

No context found.

W.W. Hager. Condition Estimates, SIAM J. Sci. Stat. Comp., 5:311-316, 1984.

No context found.

W. W. Hager, Condition estimates, SIAM J. Sci. Stat. Comput., 5 (1984), pp. 311--316.

No context found.

W. W. Hager, (1984), Condition estimates, SIAM J. Sci. Stat. Comput., 5, 311-- 316.

Online articles have much greater impact   More about CiteSeer.IST   Add search form to your site   Submit documents   Feedback  

CiteSeer.IST - Copyright Penn State and NEC