Editing Cholesky decomposition (section)

===Block variant===
When used on indefinite matrices, the '''LDL'''* factorization is known to be unstable without careful pivoting;<ref>{{cite book|last=Nocedal|first=Jorge|title=Numerical Optimization|year=2000|publisher=Springer}}</ref> specifically, the elements of the factorization can grow arbitrarily. A possible improvement is to perform the factorization on block sub-matrices, commonly 2 × 2:<ref>{{cite journal
 | last = Fang | first = Haw-Ren
 | doi = 10.1093/imanum/drp053
 | issue = 2
 | journal = IMA Journal of Numerical Analysis
 | mr = 2813183
 | pages = 528–555
 | title = Stability analysis of block <math>LDL^T</math> factorization for symmetric indefinite matrices
 | volume = 31
 | year = 2011}}</ref>

<math display=block>\begin{align}
\mathbf{A} = \mathbf{LDL}^\mathrm{T} & =
\begin{pmatrix}
 \mathbf I & 0 & 0 \\
 \mathbf L_{21} & \mathbf I & 0 \\
 \mathbf L_{31} & \mathbf L_{32} & \mathbf I\\
\end{pmatrix}
\begin{pmatrix}
 \mathbf D_1 & 0 & 0 \\
 0 & \mathbf D_2 & 0 \\
 0 & 0 & \mathbf D_3\\
\end{pmatrix}
\begin{pmatrix}
 \mathbf I & \mathbf L_{21}^\mathrm T & \mathbf L_{31}^\mathrm T \\
 0 & \mathbf I & \mathbf L_{32}^\mathrm T \\
 0 & 0 & \mathbf I\\
\end{pmatrix} \\[8pt]
& = \begin{pmatrix}
 \mathbf D_1 &   &(\mathrm{symmetric})   \\
 \mathbf L_{21} \mathbf D_1 & \mathbf L_{21} \mathbf D_1 \mathbf L_{21}^\mathrm T + \mathbf D_2& \\
 \mathbf L_{31} \mathbf D_1 & \mathbf  L_{31} \mathbf D_{1} \mathbf L_{21}^\mathrm T + \mathbf  L_{32} \mathbf D_2 & \mathbf L_{31} \mathbf D_1 \mathbf L_{31}^\mathrm T + \mathbf L_{32} \mathbf D_2 \mathbf L_{32}^\mathrm T + \mathbf D_3
\end{pmatrix},
\end{align}
</math>

where every element in the matrices above is a square submatrix. From this, these analogous recursive relations follow:

<math display=block>\mathbf D_j = \mathbf A_{jj} - \sum_{k=1}^{j-1} \mathbf L_{jk} \mathbf D_k \mathbf L_{jk}^\mathrm T,</math>
<math display=block>\mathbf L_{ij} = \left(\mathbf A_{ij} - \sum_{k=1}^{j-1} \mathbf L_{ik} \mathbf D_k \mathbf L_{jk}^\mathrm T\right) \mathbf D_j^{-1}.</math>

This involves matrix products and explicit inversion, thus limiting the practical block size.