Page 147 - Compact Numerical Methods For Computers

136 Compact numerical methods for computers
However, it is possible to use the decomposition (11.2) in a simpler fashion since
2
T
AX = ZD Z XE (11.8)
so that
- 1 T
T
T
-l
(D Z AZD )(DZ X)=(DZ X)E (11.9)
or
A Y = YE (11.10)
2
where
T
Y=DZ X (11.11a)
and
- 1
A = D Z AZD . (11.11b)
- 1 T
2
Another approach is to apply the Choleski decomposition (algorithm 7) to B so
that
T
AX = LL XE (11.12)
where L is lower-triangular. Thus we have
- 1
- T
T
T
(L AL )(L X)=(L X)E (11.13)
or
A W = WE. (11.14)
3
Note that A can be formed by solving the sets of equations
3
LG = A (11.15)
and
T
A L = G (11.16)
3
or
(11.17)
so that only the forward-substitutions are needed from the Choleski back-solution
algorithm 8. Also. the eigenvector transformation can be accomplished by solving
T
L X = W (11.18)
requiring only back-substitutions from this same algorithm.
While the variants of the Choleski decomposition method are probably the
most efficient way to solve the generalised eigenproblem (2.63) in terms of the
number of arithmetic operations required, any program based on such a method
must contain two different types of algorithm, one for the decomposition and one
to solve the eigenproblem (11.13). The eigenvalue decomposition (11.2), on the
other hand, requires only a matrix eigenvalue algorithm such as the Jacobi
algorithm 14.
Here the one-sided rotation method of algorithm 13 is less useful since there is
no simple analogue of the Gerschgorin bound enabling a shifted matrix
A'= A + k B (11.19)
to be made positive definite by means of an easily calculated shift k. Furthermore,
the vectors in the Jacobi algorithm are computed as the product of individual