Page 97 - Compact Numerical Methods For Computers

86 Compact numerical methods for computers

(Again, the diagonal elements must be chosen to be positive in the decomposi-
tion.) Equations (7.17) and (7.18) give bounds to the size of the subdiagonal
elements of L, which suggests the algorithm is stable. A much more complete
analysis which confirms this conjecture is given by Wilkinson (1961) who shows
T
the matrix LL as computed is always close in some norm to A.
Once the Choleski decomposition has been performed, linear equations
T
Ax=LL x = b (7.19)
can be solved by a combination of a forward- and a back-substitution, that is
Lv = b (7.20)
followed by
T
Rx = L x = v (7.21)
T
where we have used R to emphasise the fact that L is upper-triangular. In a
computer program, b, v and x can all occupy the same storage vector, that is, v
overwrites b, and x overwrites v. The solution of (7.20) is termed forward-
substitution because the triangular structure defines the elements v j in the order
1, 2, . . . , n, that is
v = b /L 1 1 (7.22)
l
l
and
for j = 2, 3, . . . , n. (7.23)

Likewise, the solution elements x of (7.21) are obtained in the backward order n,
j
(n – 1), . . . , 1 from
(7.24)
x = v /L n n
n
n
(7.25)

(7.26)

7.2. EXTENSION OF THE CHOLESKI DECOMPOSITION TO
NON-NEGATIVE DEFINITE MATRICES

When A is non-negative definite, the inequality (7.9) becomes
T
x Ax > 0 for all x 0 (7.27)
and inequalities (7.10), (7.12), (7.17) and (7.18) must be amended similarly.
There is no difficulty in performing the decomposition unless a zero diagonal
element appears in L before the decomposition is complete. For L mm = 0, equa-
tions (7.3) and (7.8) are satisfied by any values of L for i > m. However, if we
im
desire L to be non-negative definite and wish to satisfy the amended form of
(7.18), that is
(7.28)