If the LU decomposition of the n × n symmetric matrix A exists and the factorization of a related matrix
| $$\mathbf{A^{\prime}=A+\Delta A}$$ |
is desired, factor update is often the procedure of choice.
Section 6 examines factor update techniques for dense, asymmetric matrices. The current section examines techniques that exploit computational efficiencies introduced by symmetry. Symmetry reduces the work required to update the factorization of A by half, just as it reduces the work required to decompose A in the first place.
More specifically, the current section examines procedures for updating the factors of A following a symmetric rank one modification
| $$\mathbf{A^{\prime}=A}+\alpha\mathbf{y{y}^{T}}$$ |
where $\alpha$ is a scalar and y is an n vector.
Algorithm C1 of
Gill
| for $j=1,\cdots,n$ |
| $\delta=d_{j}$ |
| $p=y_{j}$ |
| $d_{j}=d_{j}+\alpha p^{2}$ |
| $\beta=\displaystyle\frac{\alpha p}{d_{j}}$ |
| $\alpha=\displaystyle\frac{\alpha\delta}{d_{j}}$ |
| for $i=j+1,\cdots,n$ |
| $y_{i}=y_{i}-pl_{ij}$ |
| $l_{ij}=l_{ij}+\beta y_{i}$ |
See Section 6.1 for a discussion of updating asymmetric LDU factorizations following rank one changes to a matrix.
If the LU decomposition of the symmetric matrix A exists and U is stored explicitly, the recurrence relations of the inner loop of Algorithm 20 must change. Following a train of thought similar to the derivation of Algorithm 9 (see Section 6.2 for details) results in Algorithm 21 which updates U based on a rank one change to A.
| for $i=1,\cdots,n$ |
| $\delta=u_{ii}$ |
| $p=y_{i}$ |
| $u_{ii}=u_{ii}+\alpha p^{2}$ |
| $\beta=\alpha p$ |
| $p=\displaystyle\frac{p}{\delta}$ |
| $\delta=\displaystyle\frac{u_{ii}}{\delta}$ |
| $\alpha=\displaystyle\frac{\alpha}{\delta}$ |
| for $j=i+1,\cdots,n$ |
| $y_{j}=y_{j}-pu_{ij}$ |
| $u_{ij}=\delta u_{ij}+\beta y_{y}$ |
If U is maintained in a zero-based linear array, Algorithm 21 changes in the normal manner, that is
The double subscript notation is replaced by the indexing rule defined in Equation 87.
The outer loop counter i ranges from zero to n-1.
The inner loop counter j ranges from i+1 to n-1.
The outer loops of the symmetric factor update algorithms do not have to begin at one unless y is full. If U has any leading zeros, the initial value of i (or j in Algorithm 20) should be the index of yi the first nonzero element of y. If there is no a priori information about the structure of y but there is a high probability of leading zeros, testing yi for zero at the beginning of the loop might save a lot of work. However, you must remember to suspend the test as soon as the first nonzero value of yi is encountered.
For a fuller discussion of the derivation and implementation of LU factor update, see Section 6.2.