Page 138 - Compact Numerical Methods For Computers
P. 138
Real symmetric matrices 127
that
(10.41)
(k+1) ( k )
so that each rotation causes A to be ‘more diagonal’ than A .
Specification of the order in which the off-diagonal elements are to be made
zero defines a particular variant of the Jacobi algorithm. One of the simplest is to
take the elements in row-wise fashion: (1, 2), (1, 3), . . . , (1, n), (2, 3), (2, 4), . . . ,
(2, n), . . . , (n – 1, n ). Such cyclic Jacobi algorithms have not been proved to
converge in general, except in the case where the angle of rotation f is
constrained so that
- p/4< f < p/4. (10.42)
Similarly to the orthogonalisation algorithm of §3.3, the difficulty lies in
demonstrating that the ordering of the diagonal elements is stable. For f
restricted as in (10.42), Forsythe and Henrici (1960) have carried out the
complicated proof, and most authors (Wilkinson 1965, Rutishauser 1966,
Schwarz et al 1973, Ralston 1965) discuss algorithms using angle calculations
which satisfy the restriction (10.42). In fact, among the variety of texts available
on numerical linear algebra, the author has noted only one (Fröberg 1965) which
uses the calculation based on equations (10.38), (10.39) and (3.22)-(3.27). The
advantage of the calculation given here is that the eigenvalues, which are
produced as the diagonal elements of the matrix that is the limit of the sequence
(k)
A , are ordered from most positive to most negative. Most applications which
require eigensolutions also require them ordered in some way and the ordering
that the current process yields is the one I have almost invariably been required to
produce. No extra program code is therefore needed to order the eigenvalues and
eigenvectors and there may be some savings in execution time as well.
We have already noted in chapter 3 the research of Booker (1985) relating to the
convergence of the cyclic Jacobi method with the ordering angle calculation.
Jacobi (1846) did not use a cyclic pattern, but chose to zero the largest
off-diagonal element in the current matrix. This process has generally been
considered inappropriate for automatic computation due to the time required for
the search before each plane rotation. However, for comparison purposes I
modified a BASIC precursor of algorithm 14. The changes made were only as many as
were needed to obtain the zeroing of the largest off-diagonal element in the
present matrix, and no attempt was made to remove computations present in
the algorithm which were used to provide convergence test information. The
processor time required for the order-5 Moler matrix (appendix 1) on a Data
General NOVA operating in six hexadecimal digit arithmetic was 15·7 seconds
for algorithm 14 while the Jacobi method required 13·7 seconds. Indeed the latter
required only 30 rotations while algorithm 14 took 5 sweeps (up to 50 rotations).
The comparison of timings on this one example may be misleading, especially as
the system uses an interpreter, that is, programs are executed by translating the
BASIC at run time instead of compiling it first. (This has an analogy in the
translation of a speech either in total or simultaneously.) However, for matrices
with only a few large off-diagonal elements, the relatively simple changes to