Page 61 - Compact Numerical Methods For Computers

50 Compact numerical methods for computers

This is a simpler angle calculation than that of §3.3 for the orthogonalisation
process, since it involves only one square root per rotation instead of two. That is,
if
(4.5)
then we have

c = z /p (4.6)
1
and
s = y /p. (4.7)
1
It is possible, in fact, to perform such transformations with no square roots at
all (Gentleman 1973, Hammarling 1974, Golub and Van Loan 1983) but no way has
so far come to light for incorporating similar ideas into the orthogonalising rotation
of §3.3. Also, it now appears that the extra overhead required in avoiding the square
root offsets the expected gain in efficiency, and early reports of gains in speed now
appear to be due principally to better coding practices in the square-root-free
programs compared to their conventional counterparts.
The Givens’ transformations are assembled in algorithm 3 to triangularise a real
m by n matrix A. Note that the ordering of the rotations is crucial, since an
element set to zero by one rotation must not be made non-zero by another.
Several orderings are possible; algorithm 3 acts column by column, so that
rotations placing zeros in column k act on zeros in columns 1, 2, . . . , (k - 1) and
leave these elements unchanged. Algorithm 3 leaves the matrix A triangular, that
is
A[i,j] = 0 for i > j (4.8)

which will be denoted R. The matrix Q contains the transformations, so that the
original m by n matrix is
A = QR. (4.9)
In words, this procedure simply zeros the last (m – 1) elements of column 1,
then the last (m – 2) elements of column 2, . . . , and finally the last ( m – n )
elements of column n.
Since the objective in considering the Givens’ reduction was to avoid storing a
large matrix, it may seem like a step backwards to discuss an algorithm which
introduces an m by m matrix Q. However, this matrix is not needed for the
T
solution of least-squares problems except in its product Q b with the right-hand
side vector b. Furthermore, the ordering of the rotations can be altered so that
they act on one row at a time, requiring only storage for this one row and for the
resulting triangular n by n matrix which will again be denoted R, that is

(4.10)

In the context of this decomposition, the normal equations (2.22) become

(4.11)