@prefix vivo: . @prefix edm: . @prefix ns0: . @prefix dcterms: . @prefix skos: . vivo:departmentOrSchool "Science, Faculty of"@en, "Mathematics, Department of"@en ; edm:dataProvider "DSpace"@en ; ns0:degreeCampus "UBCV"@en ; dcterms:creator "Yamamura, Eddie Akira"@en ; dcterms:issued "2011-11-21T18:57:46Z"@en, "1962"@en ; vivo:relatedDegree "Master of Arts - MA"@en ; ns0:degreeGrantor "University of British Columbia"@en ; dcterms:description """The purpose of this thesis is to give a survey of the methods currently used to solve the numerical eigenvalue problem for real symmetric matrices. On the basis of the advantages and disadvantages inherent in the various methods, it is concluded that Householder's method is the best. Since the methods of Givens, Lanczos, and Householder use the Sturm sequence bisection algorithm as the final stage, a complete theoretical discussion of this process is included. Error bounds from a floating point error analysis (due to Ortega), for the Householder reduction are given. In addition, there is a complete error analysis for the bisection process."""@en ; edm:aggregatedCHO "https://circle.library.ubc.ca/rest/handle/2429/39176?expand=metadata"@en ; skos:note "In presenting t h i s t h e s i s i n p a r t i a l f u l f i l m e n t o f the requirements f o r an advanced degree a t the U n i v e r s i t y of B r i t i s h Columbia, I agree t h a t the L i b r a r y s h a l l make i t f r e e l y a v a i l a b l e f o r reference and study. I f u r t h e r agree that permission f o r e xtensive copying of t h i s t h e s i s f o r s c h o l a r l y purposes may be granted by the Head of my Department o r by h i s r e p r e s e n t a t i v e s . I t i s understood t h a t copying or p u b l i c a t i o n of t h i s t h e s i s f o r f i n a n c i a l g a i n s h a l l not be allowed without my w r i t t e n permission. Department of Mathematics The U n i v e r s i t y of B r i t i s h Columbia, Vancouver 3, Canada. Date September 8. 1962 i ABSTRACT The purpose of t h i s thesis i s to give a survey of the methods currently used to solve the numerical eigenvalue problem f o r r e a l symmetric matrices. On the basis of the advantages and disadvantages inherent i n the various methods, i t i s concluded that Householder's method i s the best. Since the methods of Givens, Lanczos, and Householder use the Sturm sequence bi s e c t i o n algorithm as the f i n a l stage, a complete th e o r e t i c a l discussion of t h i s process i s included. Error bounds from a f l o a t i n g point error analysis (due to Ortega), for the Householder reduction are given. Iii addition, there i s a complete error analysis f o r the bi s e c t i o n process. I hereby c e r t i f y that t h i s abstract i s satisfactory. ACKNOWLEDGEMENTS The author wishes to express his gratitude to his supervisor Dr. T.E. Hu l l f o r his patient and helpfu l guidance i n preparing t h i s thesis. To the s t a f f of the Computing Centre of the University of B r i t i s h Columbia, the author expresses his thanks for the hours they spent i n helping i n the w r i t i n g of a program. F i n a l l y , the author thanks the National Research Council of Canada for t h e i r f i n a n c i a l support. i i i TABLE OF CONTENTS CHAPTER I Introduction and Summary ' 1 CHAPTER I I Methods f o r Real Symmetric Matrices 3 C l a s s i c a l Methods 3 Recent Methods k Jacobi's Method 5 Givens 1 Method 6 Lanczos 1 Method 7 Householder's Method 10 Bisection Method f o r the Eigenvalues of a Real Symmetric Tri p l e Diagonal Matrix l 6 CHAPTER I I I Error Analysis 23 Preliminaries and Notation 23 Error Bounds for the Householder Algorithm 28 Error Bounds f o r the Bisection Process 32 Conclusion 3U BIBLIOGRAPHY -36 1 CHAPTER I I n t r o d u c t i o n and Summary According t o Lanczos [ 9] , m a t r i x theory has i t s o r i g i n i n the s o l u t i o n of simultaneous l i n e a r a l g e b r a i c equations. Once a complete s y m b o l i z a t i o n of alg e b r a was introduced,, a general s o l u t i o n of a system of equations by Cramer's r u l e was discovered. However, the emphasis was s t i l l on a r i t h m e t i c . During the nineteenth century, i n t e r e s t i n the o p e r a t i o n a l aspects o f mathematics came i n t o focus. Cayley (1859) extended the realm of al g e b r a by showing t h a t a m a t r i x can be regarded as one s i n g l e a l g e b r a i c operator. The theory of the c h a r a c t e r i s t i c equation was developed by S y l v e s t e r and Wei e r s t r a s s and f i n a l l y , Frobenius gave a complete a l g e b r a i c theory. Fredholm (19OO) extended the a l g e b r a i c theory of the c h a r a c t e r i s t i c equation t o the case of i n f i n i t e l y many v a r i a b l e s , thus l a y i n g the foundation f o r the geometric treatment of l i n e a r d i f f e r e n t i a l and i n t e g r a l operators. The c h a r a c t e r i s t i c equation w i t h the a s s o c i a t e d eigenvalues and eigenvectors has many f i e l d s of a p p l i c a t i o n . These i n c l u d e v i b r a t i o n s , atomic and molecular o s c i l l a t i o n s of p a r t i c l e s , boundary value problems, and f a c t o r a n a l y s i s . E v i d e n t l y then, a knowledge of the methods a v a i l a b l e f o r the numerical s o l u t i o n of the eigenvalue problem i s important. The purpose of t h i s t h e s i s i s to give an e x p o s i t i o n of these methods f o r r e a l symmetric matrices. The essay has two main s e c t i o n s . We begin Chapter'II by d i s c u s s i n g , b r i e f l y , the determinant and \" s e r i a l \" methods f o r o b t a i n i n g eigenvalues. The shortcomings of these methods are p o i n t e d out. Then, the more s u c c e s s f u l methods of J a c o b i , Givens, and Lanczos are described i n some d e t a i l , and, we complete the d e s c r i p t i o n s by g i v i n g a d e t a i l e d account of Householder's r e d u c t i o n algorithm. Reference to d e t a i l e d accounts, p r o o f s of convergence, and e r r o r analyses are provided where a v a i l a b l e . The l a s t 2. section of Chapter I I deals with the Sturm sequence algorithm which i s used as the f i n a l stage i n the methods of Givens, Lanczos, and Householder. O r i g i n a l l y , we had planned to obtain a f l o a t i n g point error analysis for the Householder reduction and to present the d e t a i l s i n Chapter I I I . In addition, several numerical experiments were planned. Before t h i s work was completed, James Ortega's paper appeared and we discovered that he had treated, i n d e t a i l , a l l that we had planned. As a r e s u l t , i n Chapter I I I , we f i r s t present the basic preliminaries necessary for any f l o a t i n g point error analysis and then l i m i t ourselves to stating the results obtained by Ortega. For completeness of treatment, we also give an error analysis f o r the Sturm sequence bi s e c t i o n algorithm. In the l a s t section of t h i s essay, we j u s t i f y our emphasis on the Householder algorithm and indicate an area for further research. 3-CHAPTER I I METHODS FOR REAL SYMMETRIC MATRICES C l a s s i c a l Methods Here, as i n the rest of t h i s chapter, we l e t A = (a^-j), with a-y = for i , j = 1,2,...,N, denote a r e a l symmetric matrix. From a t h e o r e t i c a l point of view, i t i s apparent that the eigenvalues of A could be determined by finding the N r e a l roots of the equation det (A - X I) = 0 which i s an N-th degree polynomial equation i n X . Unfortunately,. a satisfactory r e a l i z a t i o n ( i . e . on an automatic computer) of t h i s theory i s not yet feasible f o r matrices of r e l a t i v e l y large order, say N'-IOO. The d e t a i l s i n support of t h i s statement are to be found i n a paper by H.H. Goldstine, F.J. Murray, and J. von Neumann [Id. We sketch, b r i e f l y , some of t h e i r r e s u l t s . The direct use of det (A - Kl) involves two problems. These are to determine the c o e f f i c i e n t s C^ ( i = 1,2,...,N) of the equation. det (A - X I) = A N + C X X N-X+ C 2 \\ W\" 2 +... +C N = 0 and then, to determine the N r e a l roots. Goldstine, Murray, and von Neumann divide the known methods for determining the C^ into three classes tk; p.6o] . One of these classes i s rejected on the grounds that the number of multiplications required i s of the order of 2^ [h; p.6o3 which i s a p r o h i b i t i v e figure for N~100. The methods i n the other two classes are rejected by giving an example where the r a t i o of the largest to the smallest c o e f f i c i e n t i s of the order 10^3 [h; p.6l] . They conclude [h; p.6l] : It i s very d i f f i c u l t f o r us to see how any procedure which gets a l l the c o e f f i c i e n t s C-^,. . . at one time, .'. ., can give results with any acceptable precision k. unless a very large number of d i g i t s are carried throughout. The authors next discuss the problems inherent i n root finding algorithms. Again i t i s shown C k; P-62J. that accuracy would be obtained only i f a large number of d i g i t s are carried throughout. The authors go on to discuss the \" s e r i a l \" methods f o r finding the characteristic values of a matrix. They point out that most of these methods depend upon the spectral decomposition of the matrix A. An example i s the power method-(see e.g. \\_l6^ , P«33')? The authors also point out that these methods are co s t l y i n matrix mu l t i p l i c a t i o n s and that to get an -5 accuracy of 10 i n the largest determined eigenvalue, 15 decimal d i g i t s must be carried to allow for loss of precision due to the inherent i n s t a b i l i t i e s of these methods \\_k; p.65] • Another d i f f i c u l t y of most of these schemes i s that, i n case a l l eigenvalues are desired, we must be aware of the fact that the approximation, A ^ , to the i - t h eigenvalue i s contaminated by the errors i n the previous ones - namely, X-L , \\ £ \\_k; p.65] . Thus, with regard to the c l a s s i c a l and s e r i a l methods, we believe, as does Givens V-3~\\ > that Goldstine, Murray, and von Neumann have shown that these methods are unsuitable f o r use with automatic computers i f a l l eigenvalues are wanted and the matrix i s of r e l a t i v e l y high order, say N<~100. Recent Methods We now begin the description of methods currently i n use. These methods do not require computation of the c o e f f i c i e n t s of the characteristic equation. Moreover, they y i e l d the eigenvalues i n such a manner that any error i n the approximation, X ^ to the i - t h eigenvalue of A i s not contaminated by errors i n X^, A ,^ . . . , X ^ ^ . Most of the descriptions were taken from a paper by Paul A. White £13]. The J a c o b i Method In 1 8 4 6 , C G . J . J a c o b i C T3 introduced a method of reducing A to diagonal form by means of a sequence of simple orthogonal transformations known as ( r ) plane r o t a t i o n s . I f we l e t U.. be the r - t h orthogonal m a t r i x , then we ob t a i n a sequence of transformed matrices A^ r^ w i t h A^ 0^ = A and A^ r^ = ST) ^ T - l } U » where U 1 0 u ( r ) 1 1 = cos e r u j j = s i n e r = s i n «r u j i = cos e r u k k = 1 , k = i , . . . , N = 0 otherwise, the angle, © r, bein g chosen so th a t the elements i n the ( i , j ) and ( j , i ) p o s i t i o n s become zero. That i s , J a c o b i ' s method depends on the choice of U.. such t h a t A ( r ) u^-D' ^ l / A U U ) ^ ( - 1 ) ( r ) _ _ A^ - U l r J r U l r _ 1 J r _ 1 ... U l l ( j ]_ A U X l J l .... U l r _ 1 j r _ 1 U X r J r D ( r ) where D i s a diagonal m a t r i x and Un- ^ i s the transpose of U,- * . J a c o b i j - r j r * ! r J r a n n i h i l a t e d ( i . e . r o t a t e d to zero) the maximum o f f - d i a g o n a l element at each stage, and he was able t o prove t h a t a f t e r some f i n i t e number of steps, M, a l l o f f - d i a g o n a l elements would be l e s s i n magnitude than any preassigned t > 0 . J.H. W i l k i n s o n [17] gives the d e t a i l s f o r t h i s method and a l s o discusses the p r a c t i c a l d e t a i l s f o r a c t u a l numerical work. G o l d s t i n e , Murray, and von Neumann give a thorough t h e o r e t i c a l d i s c u s s i o n of the J a c o b i method. There are two evident drawbacks t o t h i s method. The f i r s t i s t h a t scanning the m a t r i x f o r the l a r g e s t o f f - d i a g o n a l element at each stage may be time consuming. The second i s t h a t the nature of the orthogonal 6. transformations i n no way guarantees that an element once rotated to zero w i l l remain zero throughout successive stages. Hence, the scanning must be done over the entire set of off-diagonal elements. There i s , however, a complete error analysis f o r t h i s method Because of the drawbacks, White notes that t h i s method has been replaced, i n practice, by two variations which we now describe. The f i r s t of these i s known as the c y c l i c Jacobi method. This method removes the f i r s t drawback by systematically reducing to zero i n turn each element of the f i r s t row, regardless of size , provided of course, the element i s not already small enough; then, the second row, and then, the t h i r d , etc. Because of the second drawback, t h i s procedure i s it e r a t e d u n t i l a l l off-diagonal elements are s u f f i c i e n t l y small. With s u f f i c i e n t r e s t r i c t i o n s on the rotation angles, G. Eorsythe and P. Henrici C2^ have been able to prove that t h i s method converges. The second v a r i a t i o n i s r e a l l y a modification of the preceding method and i s due to Pope and Tompkins [ l 2 ^ • We start with some threshold value l > t > 0 and reduce to zero f i r s t , only those off-diagonal elements whose magnitude exceeds t. Iteration i s done u n t i l a l l off-diagonal elements are o £ ^ - Y b - . . . - ( J b — r r — r r r r+1 \" r —r-2 r — 1 such that b i s orthogonal to each of b , b ..... b -r+l & -r' - r - 1 - 1 ' This gives oC = ° Ab r ^ — r b T b — r — r b * b r = — r — r T b nb , with the other constants being zero. The l a s t two constants, Q( and p > are obtained from b = Ab ~ <-,/b - f i b N ' N - T I + 1 ~N N-N \"N - N - 1 by choosing b orthogonal to b and b f ., , —N+l ° —N — N - l . I t can be shown |_17; p . l39J that b ^ + 1 i s necessarily orthogonal to b_^ ^,. . ., b^ and that b ^ + i i s the n u l l vector. The process terminates at t h i s stage. The above description t a c i t l y assumes that no b i s the n u l l vector. — r Such an assumption i s not v a l i d \\l7; pp . l39 -l^o]] • I n case b^ i s the n u l l vector, we replace b by an arbitrary vector £^ which i s orthogonal to b , b . b } and continue the process. The only change introduced i s \" \" r - 1 ' — r - 2 ' ' — 1 that |2 =0 ; p-l^ o j . We now form the matrix B defined as B = j * 2 : •••• : \\ > or B = (5.1 I \\ \\ \\ l r \\ -••\\\\) ~ i'e-» t h e matrix B has column vectors equal to the above constructed orthogonal set of vectors. I t can be shown [ i f ; p p • I U I - I U 2 I that 9-B\" x AB = 3 3 N-1 1 P N °^N or B AB = £ 2 0 ^ Pi N-1 1 N N i f , f o r example, b^ = 0. Thus B AB i s si m i l a r to A. The eigenvalues of A are found by determining the eigenvalues of B 1 AB . We can, with a l i t t l e work, consider a symmetric triple-diagonal matrix C instead of B 1AB which i s not symmetric. This matrix C i s obtained by normalizing the vectors b. ^ 1 7 ; PP-1^7-15l\"l • Consequently, the Sturm sequence method f o r symmetric triple-diagonal matrices • also applies here. This completes our descriptions (other than Householder's method) of the methods used f o r finding the eigenvalues of a r e a l symmetric matrix. We point out that the above are descriptions only and that f o r numerical work these descriptions hardly suff i c e . We now consider i n d e t a i l the method 10-due to Householder. Householder's Method [_15~\\ Householder \\_€>~\\ suggested that the orthogonal s i m i l a r i t y transformation, used i n reducing a symmetric matrix A to triple-diagonal form, be obtained as a product of simple orthogonal matrices, P, given by the form (1) P = I = 2 w w T where w i s a column vector such that T (2) w w = 1. I t i s easy to show that P i s symmetric and orthogonal. The symmetry i s obvious, and the orthogonality then follows since (3) P TP = (I - 2 w w T ) (I - 2 w w 1 T rp = I - i + W W + ^ v W = I. In order to make Householder's method e x p l i c i t , we begin by defining a column vector w bv — r J so that w_^ i s a vector with i t s f i r s t (r-~l) components equal to zero. We then take P^ to be a P matrix as given by equation ( l ) with w = w^ . The transformation of a given ( N X \" N ) , r e a l symmetric matrix, A = (a,. ,)} to triple-diagonal form i s effected by (W-2) successive s i m i l a r i t y transformations ? 2 ' Ty-> P N - 1 • I f w e l e t A = A ^ > t h e n A ^ i s defined by the equation (5) A ( r > c P A ^ \" 1 ) P r r ( r - l ) where AK contains ( N - r ) elements i n row ( r - l ) each of which i s to be reduced to zer6 by the transformation with matrix P . This gives us (N-r) r equations to be s a t i s f i e d by the ( N - r + l ) elements of w^ . From equation ( 2 ) , 11. as applied to v_r, we obtain 2 2 Xr+1 + ••* + % 1 . These (N-r+l) equations determine the (N-r+l) elements of w_r but, because, as w i l l be shown presently, there i s a square involved, we are able to choose a determination which w i l l give the greatest numerical s t a b i l i t y or convenience. Before we consider i n d e t a i l the algebra that i s involved, we prove two simple facts about the transformation with matrix Pr. The f i r s t of these i s Result 1: The transformation with matrix P r leaves undisturbed the zeros i n rows and columns 1,2,...,r-2. This res u l t i s routine once the form of matrix P r i s shown. Evidently, r - l \\ P r -0 0 0 1 0 0 1-2X ~2X2^T+2_ -2X r + 1X r 1-2X r+1 -2XrXjj -2X r+lX N \" 2 % X r ( l \" 2 X N _ 1 ) - 2 X N _ 1 X N 2 1-2X. Vr-3 Thus premultiplying any (MSN) matrix B by P r leaves the f i r s t r-2 rows and columns of B unaltered. Wow suppose B = *1 ?2 fi2 *2 ''Pr-2 *r-2 Pr-1 0 Pr-1 X X 0 X X 12. and we post mu l t i p l y by P i t i s c l ea r that row and column ( r- l ) are the f i r s t to be a l t e r ed . Th is v e r i f i e s the f i r s t r e s u l t . The second simple f a c t i s : Resul t 2. I f A ^ = AC1-\"1) P , where A^1\"-1) and P are as be fo re , then r r> r the sum of the squares of the elements of row ( r- l ) of A^ R~^ remains i nva r i an t . Proof . Because of the above d i s cuss i on , . i t i s s u f f i c i e n t to show the r e s u l t f o r a l l a12* \" -a N l ^2 a. IN , with a . . = a . . , and P2 = l - 2 X r 0 - 2 X N X 2 \" 2 X2 X N 1-2X, •N Since P 2 A P £= (I - 2 w £ w 2 T ) A (l-2 w 2 ^ T ) = A - 2 w 2 w 2T A - (wgT A 2 ) wg T - 2 A w 2 - w 2 (w 2 T A w 2 ) ' w 2 T , (6) P 2 AP 2 = A - 2 w 2 q T - 2 q i w £ T where T ( 7 ) q = A w g - (vg A w 2 ) w2 = p - kw. with T T (8) k = (w A w ) = w p _2 _ 2 - -13-Consequently, the f i r s t row of V^kV^ ±s ll3 ~ 2 q l X3 ( a l l ' a12 _ 2 q l X2 ' a l 3 \" 2 q l X3 > •'• > a l N \" 2 q l *N } Thus the sum of the squares i s ( a ] L 1 ) 2 + ( a 1 2 - 2 q x X 2 ) 2 + ... + ( a 1 N - 2 X N ) 2 = a 1 1 + a 1 2 + ... + a 1 N t o p p p qj X| + q£ .+ + o 2 Y 2 + q x X N. 2 [a12 *1 X 2 + a 1 3 q x X 3 + . + a. IN *1 *N 2 2 - a ^ + a 1 2 + . . . + a 1 N since q.-[_ = a-j_2X2 + a-j_2 X^ + . . . + a 1 N X N and X 2 + + ... + x N = 1. This establishes the second res u l t . ( i t should be noted that Wilkinson |^ 15; p-2^J states that the sum of the squares of the elements i n any row must be invariant.) Because of Result 1 above, the d e t a i l s of the transformation with matrix P r w i l l be i l l u s t r a t e d for any stage r, i f we provide the d e t a i l s f o r = A and P r = P 2 - i.e.., the f i r s t stage. Let A = ( a ^ ) i , j = 1,2,..., N and r \"2 T _ T T l e t w 2 = w1 = (0, X 2 , X ^ , . . • •, X J J ) , so that 2 2 Xg + x 3 + + x N = 1 We wish to determine P 2 such that P 2AP 2 has zeros i n positions (l,3) > (l,N) and in. (3,l), , (N,l) . Since l e f t m u l t i p l i c a t i o n of any (N X N) square matrix by P 2 leaves the f i r s t row unaltered, P 2 AP 2 has the desired zeros i f and only i f AP 2 has the desired zeros. Thus w must be chosen accordingly. Since AP 2 = A ( l - 2 Wg WgT) = A - 2 p 2 w2T , the following set of equations must be s a t i s f i e d Ik. (9) C a 1 3 - 2 P l X 3 = 0 alk \" 2 p l Xk = 0 a IN Moreover, from Result 2 (10) s a 12 \" 2 P l X 2 = - 3 where 2 a 12 + a 13 + \"* *\" + a l N I f we mu l t i p l y equation (10) by X g and the successive equations of (9) by X 3 ' x]+> ''' ' X N r e spec t i v e l y , and then add the r e s u l t i n g equations we obta in by us ing equations(3) and ( l ) 1. (11) P l = + x 2 S 2 • Subs t i t u t i ng ( l l ) i n to (10) and so lv ing f o r X^ 2 , we have - —. (12) X = \\ 1 + f i g \\ , and, 5 \" pu t t i ng ( l l ) in to equations (9) we have (13) X k _ a l k 2 X 2 S 2 where k = 3> k, . . . , N. The upper and lower s igns go together i n equations (10), (12), and (13). From equation (10), we see that ^$2 = a 12 \" 2 pl X 2 = + (S^) ^ (lk) oi± = a 1 ] L 15-where we denote the f i n a l r esulting triple-diagonal matrix by * 1 h 0 p2 U2 /?3 h *3 Pk P N - I ° ^ N - I ^ N The above choice of signs means the X 's are not uniquely defined (we referred to t h i s before) and consequently, f o r p r a c t i c a l work, we are free to choose that sign which gives greater numerical s t a b i l i t y . Let us digress a moment to ascertain what aspects must be considered so that we obtain accurate results. We refer s p e c i f i c a l l y to a paper by CT. Fike • In t h i s paper, Fike defines the P- condition, P (A), f o r any r e a l , N-square matrix A and i t s proper value o(^ • He says that \" P^ (A) can be regarded as a measure of the p r a c t i c a l d i f f i c u l t y attached to the problem of computing the proper value .«« Using t h i s P- condition, Fike goes on to show that r e a l symmetric matrices are w e l l conditioned - i . e . , there i s not too much d i f f i c u l t y i n computing a proper value o(, of A and that \" s i m i l a r i t y transformations K. made with orthogonal matrices cannot cause a deterioration i n the conditioning of the problem.\" Fike also refere to Householder and Householder and Bauer [^6j who suggested that orthogonal s i m i l a r i t y transformations are p a r t i c u l a r l y stable i n numerical work. Wilkinson [\"15] also says that i t i s essential that the matrices, P , be as accurately orthogonal as possible. As Wilkinson goes on to point out, t h i s means that, since we determine X X^ ,. . •, X^ , for example, by dividing by Xg, we should choose the sign 16. 2 i n equation (12) such that Xg i s a s large as possible. I f we do t h i s , the resu l t i n g equations are (15) (16) 1 + a12 S G N ( a l 2 ^ T\" S 2 ] Xk = a l k SGN(a l g) 2 X S 2 2 with Xg = \\j Xg2 ) since the sign i s not important, and (17) j 3 2 = - SGN(a ) where SGN(a 1 2) = f +1 i f a 1 2 > 0 -1 i f a 1 2 < 0 . I f we use equations (12) and (13) then, according to Wilkinson ^ 15; P-2U^ , the equation 2 2 2 x 2 . + X + ... + X n = 1 i s very accurately s a t i s f i e d . As a f i n a l p r a c t i c a l ' d e t a i l , we point out that the (2) (l) transformed matrices A^ , A w , (6), (7), and (8). , A (N-l) are obtained by using equations Bisection method f o r the Eigenvalues of a Real Symmetric Triple-Diagonal Matrix I t i s evident from the above descriptions of the methods of Lanczos, Givens and Householder, that we need a method f o r obtaining the eigenvalues of a r e a l symmetric triple-diagonal matrix C = (c. . ) , where 1J and c. . ... i , i + l 0 | i - j . l > i : i + l , i = P i+1 i = 1,2, , N-l 17-We s h a l l consider i n d e t a i l the bisec t i o n method f o r computing the eigenvalues of the matrix C. This method depends upon the Sturm sequence associated with the matrix (C - X i ) , X r e a l . Ortega [ ] l l ] considers the theory i n d e t a i l . He even gives a fcrivial example to show that the theory i n Givens' Oak Ridge paper i s not quite correct ^H; P-26^ . We consider the special case where none of the j^ -^'s are zero. In case there are any fi> ±'s exactly zero, they are replaced, i n the program w i i t t e n , by the smallest p o s i t i v e number recognized by the machine. According to Wilkinson QL<7; P-13O] , exactly zero P-^ 's a r e very rare and i n case there i s a ^ = 0, he fe e l s that i t i s not worthwhile to separate C into two matrices. Givens [13; p.UOl] has apparently shown that such a change can cause a change i n the eigenvalues of not more than twice the magnitude of the non zero term replacing the zero. As the error analysis given below w i l l show, such an error i s indeed of no consequence to the accuracy of the method. For the matrix C - X I and i = 0, 1, 2, ..., N, we consider a sequence f^ of the upper l e f t p r i n c i p a l minors defined by: (17) f. ( X ) = 1 i f i = 0 ( o < 1 - X ) i f i = l (°o - 1 i f f N ( \\ ) < 0 i f f 1 (X ) < 0 , / « a s f [ f N - i ] « f H <*> - 0 18. D e f i n i t i o n 2. Let A (X ) denote the number of agreements i n sign of the sequence \\^f± (X )\"j ( i = 0, 1, . .. , N ) calculated by means of d e f i n i t i o n 1. Theorem: Let a r e a l symmetric triple-diagonal matrix be given by C where none of the f^j_'s are zero. Then, f o r any r e a l Xjthe number of eigenvalues of C that are greater than X i s given by A ( X ) . Proof: We f i r s t establish the following properties of the sequence | f . (X : (a) Two consecutive f ± ( X ) cannot both vanish for the same \\ .. (b) I f f. ( X ) = 0 , then f i _ 1 ( X ) • f 1 + 1 ( \\ ) < 1 f o r i = 1, ... , N - 1. (c) f ^ ( X ) = 0 has no multiple roots. Property (a) i s proved by induction. Ob viously, both TQ and f-^ cannot be zero. Hence, we assume the res u l t f o r k i s a c l a s s i c a l Sturm sequence f o r X, not a zero of f^. That t h i s i s the case -w i l l be proved l a t e r . Thus the number of zeros of f ^ greater than X i s equal to the number of variations i n sign of (*). I f we rewrite (*) as ^ f Q > - f-j_ j + fg > \" 3 ^ then there i s an agreement of signs i n | f ^ i f and only i f there i s a corresponding difference of signs i n (*). Hence, i f X i s not a zero of f j j , A ( X ) gives the number of roots of f ^ (\\) = 0 that are greater than X . Now suppose XQ i s a zero of f^. Then by property (b) and the fact that f i s a continuous function, we may conclude that there i s some £- neighborhood of X Q , say N (. X Q , f e ) , such that f ^ (X ) £ 0 f o r X £ N ( A 0 , € ). Consequently, SGN j^f N_ x (A )] i s constant f o r X £ N ( XQ, e ) . Moreover, since fjj_]_ (AQ ) /0> the number of agreements i n sign i s constant for the sequence ( A ) , ••• , f Q ^ provided \\ € N ( A Q , )• Since X Q i s a simple zero °f fjj (X ), A ( AQ - & ) - A ( A Q + S ) = 1 for a & such that 0< S < 6 . Therefore, SGN £ f N ( X 0 - S ) .= SGN J^N-l ( ^ 0 \" ^ )J • In the l i m i t , as \\})-*Q , we have, using d e f i n i t i o n r-l, that SGN £ f N ( A Q)| / SGN ( X 0 ) J • Thus A ( X Q ) = A ( A Q- o) - 1. That i s , with our choice of signs, A (X ) gives the number of eigenvalues of C greater than A . 20. To complete the proof, there remains to show that.the sequence (*) i s , for not a root of f N ( X) = 0, a Sturm sequence. This i s v e r i f i e d by-induction, but before doing so, we give an i l l u s t r a t i o n for the sequence f^0-> \" ?±> ?2> ~ ^ \" Le't u s consider the following diagram where the horizontal l i n e s represent the X -axis and the heavy v e r t i c a l bars the roots of f i ( i = 0, 1, 2, 3, k). F 0 ± -f l ^ - ^ f + + + + + | | + + + + + + . f ^ | + + + | + + + + + + + + \\ | + + + + | | + + + + Clearly, i f the zeros of f are i n the r e l a t i v e positions shown, ^ f o > - f]_> %2> ~ ^ 3' ^ s a Sturm sequence. Consequently, f o r the general sequence (*), we need only show that t h i s r e l a t i v e positioning i s necessary. For the i n i t i a l induction step we consider ^ V Q , - f - ^ - From the above diagram, t h i s i s c l e a r l y a Sturm sequence. Assume that f o r K< N, K even, <^ f Q , - f ]_, .... , - f R _ 1 , f K ^ has the r e l a t i v e positioning mentioned. We wish to show that <^0, ,• - , *K, ~ f > - fK+1^) also has the desired r e l a t i v e positioning of the roots of f ^ - 0 ( i - 1, 2, . . . , K + 1). By the induction hypothesis, the zeros of f ^ are positioned r e l a t i v e to those of -fj^_3_ as follows: 21. -f K-1 ' J + + + | I + + + + I I + + + fK + + | - - _ l + + + | - - - l + v + + | - - - | + + + + | - _ - + + + \" fK+l : • |+ + + — 1 ifcn 1 > Consider the largest root, X of - f (X ) = 0. For X large enough 1 ' K+1 - ( X) y 0. Hence, to the right of X ^ + ^ , t h e sign of - f^- +-j_ i s p ositive. For A W w e h a v e ^ [\"_ fK_± (\\W . ) \" ] • [ - f K + 1 ( X ( f )] 0, f K + 1 ( X( K ) )< 0, ... , f K + 1 (X(K)>0, That i s , there are K-1 d i s t i n c t zeros of f K + 1 i n the i n t e r v a l j ^ X ^ , X ^ J • Clearly, the only d i s t t i b u t i o n s a t i s f y i n g the conditions i s obtained i f two consecutive zeros of f straddle a zero of fv- The argument for K odd i s K+1 K obtained by an obvious v a r i a t i o n i n the above argument. This completes the induction and we may conclude that (*) i s a Sturm sequence. Moreover, the proof of the theorem i s now complete. We now give a description of the bis e c t i o n procedure. We assume that the maximum modulus of the ^ . 's and ^ ' s i s less than 1 - i . e . , max | Q( | , | |J t | 1, where P i s given by maxjjfij + W i l + | P>i+ll| ( i = 1 , 2 , . . . , N). By the above theorem,^ i s S U c h that P-k<\"A(^ P-k + 1. The-interval of unit length (P-k, P-k + l ] i s divided into two by adding \\ to P-k. Let us put P-k + 1 = \\ so that the i n t e r v a l we consider i s («\\ -1, Xj . To determine whether Xj £ (X -1, X - ^ or (X - \\, X1| } we evaluate A(X - •§•)• Suppose X, £ (A - \\, X j — i.e., A (X - |r)>l. We now divide (X- i , XJ \"by adding (^)2= ^ to A - \\ and evaluate A (\\ - ) to determine whether A ( 6- (\\ - \\, X - 3_ J or X, ^ (X - ^ , Xj • Continuing t h i s process, we see that at the j - t h stage X^ ( l ) < A, £ A^(l) where A^Cl) - X^(l) = {^)3. In the program written, we stop at a stage j=J such that (ir) J ^ 2 . i o - * where t i s the number of d i g i t s c a r r i e d i n the mantissa of the f l o a t i n g point number format. To continue the process, we now put A^Cl) = P and repeat the above procedure requiring, of course, that the choice of intervals be made by having A (X)^>2. Thus to straddle A^by appropriate A^ (r) and \\ j ! r ) , we begin by putting A ^ r - l ) = P and require that our choice of intervals depend upon A (A) being greater than r. The process ends when we have straddledX ^ . A word on our scaling i s i n order. I t should be clear t h a t . i f we do not scale then, instead of subtracting 1 from a P-value, we would have to subtract a 1 scaled by 1 0 s , where s i s determined so that i n P - 1 0 sl we are subtracting a 1 from the f i r s t decimal d i g i t of P. And, instead of adding powers of \\, we would be adding 10 s f o r j - 1 , 2 , . . . . Since the i n i t i a l scaling method i s simpler to code, i t ' s use was adopted. The r e l a t i v e merits of the two methods beyond the coding were not considered. This completes our descriptions of the more successful methods currently-in use. We next give the results of Ortega's error analysis of the Householder re-duction and following t h i s , we give a complete error analysis of the bisect i o n procedure. 23-CHAPTER I I I ERROR ANALYSIS Prel i m i n a r i e s and Notation The notation used follows that of Wilkinson £l^3 • Exact mathematical operations w i l l be denoted by t h e i r usual symbols: \" — \" , nx\", and \" The corresponding f l o a t i n g point operations are denoted by f l (x - y ) , f l (x + y ) ; f l ( x • y ) , and f l (x y) or f l ( x/ y) where x and y are f l o a t i n g point d i g i t a l numbers. The f l o a t i n g point arithmetic subroutines used give the r e s u l t of each operation as the c o r r e c t l y rounded standard f l o a t i n g point number. This implies the following r e l a t i o n s f o r the r e l a t i v e errors that are introduced: . . . . . _. f l (x + y) = x (1 + £ ) + y ( 1 + €) f l (x.y) = x.y (1 + L) x ( i + *) • y x (I + . y ( i + f i ( x/y) = [ x / y ] ( i + 6 ) x ( l + € ) / y = x / y ( i + ^ ) where \\&\\ ^ \\ 1 0 1 * with t being the number of d i g i t s being used i n the mantissa of the f l o a t i n g point representation. I t should be noted that separate uses of an 6 do not denote, nece s s a r i l y , the same thing. For example, the meanings of the two €•1 s i n the a d d i t i o n equation are r e l a t e d only by the requirement that both be bounded above i n magnitude by \\ 10\"*\" * . Each of the above r e l a t i o n s implies that the r e s u l t of each f l o a t i n g point arithmetic operation on two f l o a t i n g point d i g i t a l numbers x and y i s the exact r e s u l t f o r one or two s l i g h t l y modified numbers. For an extended product we have f l ( x x - x 2 x N ) = x x . x 2 x N (1-:+^) (1 + € 2 ) . . . ( l +6^ 2k. with ^ ^ 1 0 1 _ t f o r i = 1, 2, ... , N. That i s , the extended product of N f l o a t i n g point numbers i s an exact result f o r N s l i g h t l y modified numbers x i ( l + Although f o r an extended sum there are no useful bounds f o r the expression f l ( x 1 + Xg + ... + x w) - 1 X X + Xg + . . . + XJJ with the additions done naturally from l e f t to r i g h t , i t i s s t i l l true that the calculated sum i s the exact sum of modified numbers x-j_, Xg, - • • , XJJ d i f f e r i n g from the corresponding x i , Xg, ... , x w by small r e l a t i v e errors. To v e r i f y t h i s we start by assuming that f l ( x 2 + xg, + ... + N ) = X l (1 + 4 N ) ) + ...+ x„ ( 1 + = A , say. Then f l ( x x + x 2 + ..• + x N +x N + 1) = f l (A + x w + 1 ) = A (1 + e ) + x N + 1 (1 + t) and a l i t t l e computation shows that ( i - | i o 1 \" * ) 1 1 £ i + € 1 ^ N + 1 ) £ ( l + \\ lo^f and 'for i = 2, 3> • • • , N+l (1 - | io 1 - t ) N + 2 - i £ 1 + € (N+l) ^ d + i 1 0 l - t ) l N + 2 - i • That i s , f l (x, + x 2 + ... + x ^ ) = x x (1 + + x 2 ( l . + ^ g ( N + l ) ) + -(N+l) , /-(N+l)v % ( i + £ \\ U x N + 1 ( i + e l N + i ) Consequently, i f we assume that the operations take place i n the natural order, the l a s t equation implies that the calculated sum i s the exact result f o r numbers x^ such that Xj_ = x. . ( 1 + £(N+1)) i = i , 2, , N+l • 25-For an inner product c a l c u l a t e d i n single p r e c i s i o n f l o a t i n g point, we have f l ( x x . y-L + x 2 ,y 2 + + XN ' ^ N^ = x, . y, (1 +4N> ) + x 2 . y2 (1 + ) +...+ ^ • yN (1 +4N)) with , . f..\\ ( i - \\ i c ^ V A^ . A calculated A There are two conditions that must be met. F i r s t , the elements of A are to be related to the corresponding elements of A so that a\\y d i f f e r s from a^j by simple multiples of a rounding error. Second, the algorithm must not be altered. There are some d i f f i c u l t i e s and we now i l l u s t r a t e two of these. The errors introduced i n c a l c u l a t i n g x + a 1 2 SGN ( a 1 2 ) 1_ S 2 cannot be attributed to the elements a.. . of A i n a simple straightforward 1 J i manner. Let us see why. In accounting f o r the error made i n computing S 2, 27-we may say that i f we used a-^ ( l + € ) for i = 2 , 3> •••> N then ^ c a l c u l a t e d S 2 ~] i s an exact re s u l t f o r these modified numbers. Let us put.. T = ajg SGN (a-^g) • Since the \\ and 1 occurring i n the formula for [\"calculated S2~\\ 2 Xg i s part of the algorithm, we want to attribute a l l errors made i n calcu l a t i n g x| = \\ j \\ + TJ to the quantity T - i.e. calculated Xg = ^1 + TJ. The errors occur when we divide |.a-|_g,| by ^calculated S 2 J , then when t h i s result i s added to 1 and f i n a l l y , when the l a s t result i s m u l t i p l i e d by \\. The f i r s t and l a s t operations give r i s e to the following problem. I f we modify a-^ g to account for these errors then we have a 1 2 ( 1 + 6 ) € ) a 1 2 (1.+ 6) 2 +. and we want to claim that the computation i s exact for a-^ g so modified. Since the three E's may stand fo r three different nonzero quantities, some objection to the claim fo r exactness may be made. We do not know whether t h i s objection i s v a l i d or not. From the expression for Xg, i t i s clear that T ^ 1. Because of t h i s f a c t , i t i s no longer true that the error made i n adding ^calculated T^Jfto 1 may be attributed to ^ calculated TJ i n such a fashion that the error i s bound by an £ . As an example, suppose we are carrying k d i g i t s and that a 1 2 | i s .0009. Then the exact re s u l t i s 1 .0009, whereas the [^calculated S 2 J calculated r e s u l t i s 1 .001. Thus, the error i s .0001 and we want to say that .0001 = .0009 ( l + £ ), where |e| ^ 5.10\"^. In fa c t , £ = - 8/q. Because of these d i f f i c u l t i e s , the \"backwards\" error analysis was abandoned. 28. Error bounds for the Householder Algorithm The p r i n c i p l e of Ortega's analysis i s based upon the following observations of Wilkinson • Using the notation introduced i n describing Householder's method l e t A^1\"^ be the eigenvalues of [ c a l c u l a t e d A^r^j with the eigenvalues of A(!) = A. Let | A . ( r ) - ^ S I ( R ) • Then for the triple-diagonal matrix [calculated A(N-1)J , we have that ' r=I . r=l The problem now becomes one of obtaining bounds for Sj^r^ • In order to do t h i s , l e t Qr+-]_ be the exact orthogonal matrix which would be derived by the Householder algorithm applied to [ c a l c u l a t e d A^r^3 • For r = 1 2, . . ., N-2 , we define E r = [calculated ( P r + i \\_ calculated A ^ } P r+i )J - ^exact ( Q r + 1 ^calculated A ^ ] Q r + 1 ) J = ^calculated A ^ r + 1 ^ j - j\"exact ( Q r + 1 [ c a l c u l a t e d A ^ J Q r + 1 )j • Then by d e f i n i t i o n , the eigenvalues of [calculated A^r+^\"^J are ^ and, by d e f i n i t i o n and s i m i l a r i t y , the eigenvalues of exact (Qr+}_ L~ ca^- A r^^ J Qr+2_ ) \\ (r) are A . . Thus, by L i d s k i i ' s theorem ) ' i (r+1) _ ^ ( r ) | ^ ^ . max eigenvalue of E r , and now the problem i s to f i n d bounds f o r the elements of E r , r = 1,2,..., N-2-Ortega obtains bounds f o r max | A-^\"^ - \\i^ n\"''\"^ i r e l a t i v e to both the spectral norm of A and the Euclidean norm of A, where these norms are defined by: || A || g = (max A± )* with A . , being the eigenvalues of A A, and 29-The thoroughness of Ortega's analysis may be i l l u s t r a t e d by the following points. Since the norms; either spectral or Euclidean, of {^calculated A ^ ^ J , .... , ^ c a l c u l a t e d A ^ N \" ^ J do not necessarily remain that of A , Ortega makes a study of the growth of the norms of [c alculated A< r >] i n terms of the norm of A . He also carries a l l higher order terms ( i . e . terms involving at least £_2) u n t i l the f i n a l stage and then finds a bound.for them. F i n a l l y , he considers both normal and double precision f l o a t i n g point inner products f o r vectors. Before we give Ortega's r e s u l t s , a few words on his notation are i n order. He denotes by m^ , ms, and nip bounds for the r e l a t i v e errors* i n the basic arithmetic operations, square roots, and inner products, respectively. The res u l t s f o r the spectral norm are: I f mg < 211^ , Nm^ £ 10~k, N2!^ 4 10\"3, ^ IO - 2, where N i s the order of A , and £ i s the maximum error i n any eigenvalue then: 1 £\\ ^ 55(N-2)ms + (3.2 N 5/ 2 + 9-75 N 2 + 6.0 N 3/ 2 + 157-0 N-397)mb , \" A \" s ^ l - 55(N-2)ms - (3-2N5/2 + 9.75 N 2 + 6.0 N 3 / 2 + 157.ON-397H f o r normal f l o a t i n g point inner products.. I f ms ^ 21%, Nm^ $ 10'k , N 3/ 2 1% ^ I O - 2 then • 55(N-2)ms + (6.0 N 3/ 2 + l 6 l . l N - .T)mb V s • • • \" -1 • • — —• 1 - -\" A \" S ^ ~ 55(N-2)ms + (6.0 N 3/ 2 + 161.1 N - 31+8.7)1% f o r accumulated inner products. The corresponding results f o r the Euclidian norm are: Ifc 1 / 55.5(N-2)ms + (13.9 N 2 + 160.9N - 378)1% || A,j ^ 1 - 55-5(N-2)ms + (13.9 N + 160.9 N - 378)1% f o r normal f l o a t i n g point inner products; and, Ifet < 55-5(N-2)ms +17^-8 (N-2)mb for accumulated inner products. || A 1| ^ 1 - 55-5 (N-2)m + 17^ -8 (N-2): im (* From our e a r l i e r discussion we should keep.in mind that f o r addition and inner products we do not generally have true r e l a t i v e errors.) 30.-Denoting the r i g h t hand side of these bounds by F (N, m^ , m s ) > Ortega has prepared the following tables f o r comparisons. Table 1 [lO; p.38^ —12 F(N,mj3 , mg ) for m^ = 5 x 10\" and mg = 2m-j3 (Spectral Worm) N F (N, mb > ms ) Normal Accumulated 10 2 10 X 10\" -8 1.18 X 10' -8 30 1 63 X 10\" •? h.31 X 10' -8 50 k 82 X 10\" -7 7.65 X 10' -8 100 2 25 X 10\" -6 1.6k X 10' -7 200 1 13 X 10' -5 3-56 X 10' -7 500 1 03 X 10\" -k 1.03 X 10\" -6 1000 5 50 X 10' •k 2.31 X 10' -6 Table 2 [10; P-38J F(N>mb , ms ) for mb = 5-x 10\" and ms = 2mb (Spectral Norm) N F (N, , m ) s ' Normal Accululated 10 2.10 x IO\"1* 1.18 x 10~k 30 1.63 x 10\"3 k.31 x 10_1+ 50 i^ .82 x 10\"3 7.65 x 10~h 100 2.25 x IO\" 2 1.6k x 10\"3 200 * 3.56 x I O - 3 500 * 1.03 x IO\" 2 1000 * 2-35 x I O - 2 * No figures of accuracy Table 3 \\ J L O j p. 53] F(N,mb, ms ) for m^ = 5 x 10~ 1 2and mg = (Euclidean Norm) N F (Nj mb> ms ) Normal Accumulated 10 •1 76 X 10\" -8 •1.1k X 10' -8 30 1 01 X 10\" •7 k.02 X 10' -8 50 2 h3 X 10\" -7 6-90 X 10' -8 100 8 32 X 10' -7 l . k l X 10' -7 200 3 07 X 10' -6 2.8k X 10' -7 500 l 82 X 10' -5 7-15 X 10' -7 1000 7 10 X 10' -5 . l.kk X 10 -6 Another p r a c t i c a l use of the bound F(N,mb, mg) i s i l l u s t r a t e d i n the following table. Table k [_10; p.Uo] Maximum allowable N so that F(N,m , m )<^ S f o r m-fr = 5 x I O - 1 2 , ms = 2m^ if Normal Accumulated 23 6k IO\" 6 69 k90 10-5 190 3.2 x 103 10'k U90 1.8 x 10k 10\"3 1270 9.k x 10^ 32. Error bounds f o r the Bisection Process To complete the error analysis, we now obtain error bounds f o r the computed eigenvalues of the symmetric triple-daigonal matrix C. The analysis i s e s s e n t i a l l y that of Wilkinson [lh; pp.324-3263 • We r e c a l l that the elements of C were sealed so that for a l l i W-J ^ 1 and |J3jJ ^ 1; also, none of the fS-ps were zero. Referring to the description of the Sturm sequence bisect i o n method, we see that a sequence ^?Q> f-j_> • • • > f-^) i s calculated. We s h a l l show that t h i s sequence i s an exact sequence for a modified matrix ( C - \\ l ) . Hence, using L i d s k i i ' s theorem, we obtain abound for A'j_ \"X^ where X^ and \\ i are the eigenvalues of C and C respectively. Since fr(K) = (o since the re l a t i o n |X|^max |a | holds fo r any matrix A. Consequently, | « ^ „ - ^ | ^ 4 for a l l t r i a l values of X. Thus, \\o(' - t^J = |(o Wilkinson compared Givens' method with one which used elementary-s i m i l a r i t y transformations and found that on the matrices (up to order 30) that were tested, the elementary transformation method was j u s t as accurate. But he points out [l5> P - 2 6 J that the error analysis i n d i c a t e d that the u n i t a r y transformations are more stable numerically. Also, C.T. Fike's p a p e r > mentioned before, indicates that such should be the case. Consequently, since Householder's reduction r e t a i n s the mentioned advantages even on the unsymmetric case, i t may turn out to be a very important'method. Research along these l i n e s i s planned f o r the future. 36-BIBLIOGRAPHY £lj Fike, CT., Note on the p r a c t i c a l computation of proper values, J. Assoc. Comp. Mach., vol.6, 1959, pp.360-362. \\_2~\\ Forsythe, G. and H e n r i c i , P., The c y c l i c Jacobi :method f o r computing the p r i n c i p a l values of a complex matrix, Trans. Am. Math. S o c , vol.94, n o . l , Jan.i960, pp.1-23• [3~\\ Givens, W., Numerical computation of the c h a r a c t e r i s t i c values of a r e a l symmetric matrix, Oak Ridge National Laboratory Report No. I57I+, 1954. \\k~\\ Goldstine, H.H., Murray, F.J., and von Neumann, J . , The Jacobi method f o r r e a l symmetric matrices, J . Assoc. Comp. Mach., vol.6, Jan. 1959, pp.60-66. Householder, A.S., Hedrick lectur e s d e l i v e r e d at the summer meeting of the Mathematical Association of America, August, 1958, c i t e d i n Fike, C.'T. [6-J Householder, A.S.,.and Bauer, F.L., On c e r t a i n methods f o r expanding the c h a r a c t e r i s t i c polynomial - Numerische Matematik, v o l . 1, n o . l , 1959, PP-29-37• Jacobi, CG.J., Uber ein l e i c h t e s Verfahren, die i n der theorie der sakular stbrqngen vorkommenden gleichungen numerisch augzulosen, J . Reine Angew. Math., 30 (1846), pp.51-95, c i t e d i n White, P.A. 8^ Lanczos, C , An i t e r a t i o n method f o r the s o l u t i o n of the eigenvalue problem of l i n e a r d i f f e r e n t i a l and i n t e g r a l operators, J . Res. N.B.S, vol.45, 1950, pp.255-282. |9 J Lanczos, C , Applied a n a l y s i s , Prentice H a l l , New Jersey, 1956, 1 J PP.49-51-TlOj Ortega, J.M., An er r o r analysis of Householder's method f o r the symmetric eigenvalue problem, Tech. Report No.18, Appl. Math. Stat. Lab., Stanford U n i v e r s i t y , 1962. 111J Ortega, J.M., On Sturm sequences f o r t r i d i a g o n a l matrices, J . Assoc. J Comp. Mach., vol.7, J u l y i960, pp.260-263. fl2 J Pope, D. and Tompkins, C , Maximizing functions of r o t a t i o n , J . Assoc. Comp. Mach., vol.4, 1957, pp.459-466. [13 J White, P.A., The computation of eigenvalues and eigenvectors of a J matrix, J . Soc. Ind. Appl. Math., vol.6, 1958, pp.393-^03-Jl4 Wilkinson, J.H., Erpor analysis of f l o a t i n g point computation, Numerische Matematik, vol.2, no.5, i960, pp.319-340. J15 J Wilkinson, J.H., Householder's method f o r the s o l u t i o n of the algebraic eigenproblem^ The Computer Journal, vol.3, A p r i l i960, pp.23-28. 37' [j-6J Wilkinson, J.H., Notes on p r a c t i c a l methods of solv i n g l i n e a r systems and c a l c u l a t i n g the eigensystems of matrices, 'National P h y s i c a l Laboratory, Teddington, England, 1959-\\j-7~\\ Wilkinson, J.H., Theory and p r a c t i c e i n l i n e a r systems and the determination of c h a r a c t e r i s t i c values and c h a r a c t e r i s t i c vectors, Summer Session (1958) Univ. of Michigan, Ann Arbor, Mich., pp.Ill - 1 5 2 . [l8^j Wilkinson, J.H., S t a b i l i t y of the reduction of a matrix to almost t r i a n g u l a r and t r i a n g u l a r forms by elementary s i m i l a r i t y transformations, J . Assoc. Comp. Mach., v o l . 6 , 1959, PP-336-359-"@en ; edm:hasType "Thesis/Dissertation"@en ; edm:isShownAt "10.14288/1.0080575"@en ; dcterms:language "eng"@en ; ns0:degreeDiscipline "Mathematics"@en ; edm:provider "Vancouver : University of British Columbia Library"@en ; dcterms:publisher "University of British Columbia"@en ; dcterms:rights "For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use."@en ; ns0:scholarLevel "Graduate"@en ; dcterms:title "Methods for the numerical solution of the eigenvalue problem for real symetric matrices"@en ; dcterms:type "Text"@en ; ns0:identifierURI "http://hdl.handle.net/2429/39176"@en .