7.4 Doolittle’s Method for Symmetric Matrices

If A is a symmetric n × n matrix, Algorithm 10 computes – one row at a time – the upper triangular matrix U that results from a Doolittle decomposition. The upper triangular portion of A is overwritten by U.

Algorithm 10: Symmetric Doolittle’s Method
for $i=1,\cdots,n$
   for $j=1,\cdots,i-1$
      $w_{j}=\displaystyle\frac{a_{ji}}{a_{jj}}$
   for $j=i,\cdots,n$
      $\alpha=a_{ij}$
      for $p=1,\cdots,i-1$
         $\alpha=\alpha-a_{ip}w_{p}$
      $a_{ij}=\alpha$

The algorithm uses a working vector w of length n to construct the relevant portion of row i from L at each stage of the factorization. Elements from L are derived using Equation 81.

If the upper triangular portion of A is stored in a linear array, Algorithm 11 results in the same factorization (assuming zero based subscripting). The upper triangular portion of A is overwritten by U.

Algorithm 11: Symmetric Doolittle’s Method - Array Based
for $i=0,\cdots,n-1$
   for $j=0,\cdots,i-1$
      $\mathbf{w}[j]=\displaystyle\frac{\mathbf{a}{[jn-j(j+1)/2+i]}}{% \mathbf{a}{[jn-j(j+1)/2+j]}}$
   for $j=i,\cdots,n-1$
      $\alpha=\mathbf{a}{[in-i(i+1)/2+j]}$
      for $p=1,\cdots,i-1$
         $\alpha=\alpha-\mathbf{a}{[in-i(i+1)/2+p]}\cdot\mathbf{w}{[p]}$
      $\mathbf{a}{[in-i(i+1)/2+j]} = \alpha$

For the general implementation of Doolittle’s method, see Section 4.2.