Introduction

Every so often in mathematics, a structure hides a shortcut. Cholesky decomposition is one of the cleanest examples: if you know a matrix is symmetric (A=ATA = A^T) and positive-definite (all eigenvalues positive, or equivalently xTAx>0x^T A x > 0 for all x≠0x \neq 0), you can factor it as

A=LLTA = L L^T

where LL is a lower-triangular matrix. Because LL and LTL^T carry the same entries — just mirrored — you only need to compute half as many numbers as you would for a general LU decomposition.

The method is named after André-Louis Cholesky, a French military officer and geodesist who devised it around 1910 while adjusting survey measurements for the map of Crete. He never published it himself; the algorithm reached the world through a posthumous note by his colleague Commandant Benoit in 1924.

Today Cholesky decomposition is one of the workhorses of numerical linear algebra: it solves systems Ax=bAx = b efficiently, powers Kalman filters, Gaussian process regression, and the "square root" step inside countless optimization algorithms.

Build L Column by Column

Below is a 3×33 \times 3 symmetric positive-definite matrix AA. The demo computes the Cholesky factor LL one column at a time, highlighting each entry as it is calculated using the recurrence formulas. You can edit the diagonal entries to try different matrices, then verify that LLT=AL L^T = A.

<!-- Cholesky demo: show A and L side by side, step through the factorization -->
<div class="container">
  <div class="section">
    <h3>Matrix A (edit to try your own)</h3>
    <p class="hint">Click any cell to change a value. A is kept symmetric automatically.</p>
    <div class="matrix-wrap">
      <table id="matA" class="matrix"></table>
    </div>
  </div>
  <div class="section">
    <h3>Cholesky factor L</h3>
    <p class="hint">Built column by column. Yellow = current entry, green = done.</p>
    <div class="matrix-wrap">
      <table id="matL" class="matrix"></table>
    </div>
  </div>
</div>
<div class="btns">
  <button id="btnStep" type="button">Next entry</button>
  <button id="btnAll" type="button">Compute all</button>
  <button id="btnReset" type="button" class="ghost">Reset</button>
</div>
<div id="status" class="status"></div>
<div id="formula" class="formula"></div>

/* styles for the Cholesky decomposition demo */
* { box-sizing: border-box; }
body { font-family: system-ui, sans-serif; color: #222; margin: 0; padding: .4rem; }
h3 { margin: 0 0 .3rem; font-size: 1rem; color: #1d3557; }
.container { display: flex; gap: 1.5rem; flex-wrap: wrap; margin-bottom: .6rem; }
.section { flex: 1 1 140px; }
.hint { font-size: .82rem; color: #555; margin: 0 0 .4rem; line-height: 1.4; }
.matrix-wrap { display: inline-block; border-left: 2px solid #1d3557; border-right: 2px solid #1d3557;
               padding: 4px 8px; }
.matrix { border-collapse: separate; border-spacing: 4px; }
.matrix td { width: 44px; height: 36px; text-align: center; font: 600 14px ui-monospace, monospace;
             background: #e8eef3; border-radius: 6px; border: 1px solid #cdd9e3; color: #1d3557;
             transition: background .25s, color .25s; }
.matrix td.editable { background: #f0f4f8; cursor: text; }
.matrix td.zero { color: #aaa; }
.matrix td.active { background: #ffd166; border-color: #e0a800; color: #333; }
.matrix td.done { background: #c8f0d4; border-color: #5cb87a; }
.matrix td.empty { background: #f7f7f7; color: #bbb; border-color: #ddd; }
.btns { display: flex; gap: .5rem; flex-wrap: wrap; margin-bottom: .5rem; }
button { font: 600 14px system-ui; padding: .4rem .85rem; border: 1px solid #1d3557;
         background: #1d3557; color: #fff; border-radius: 8px; cursor: pointer; }
button:disabled { opacity: .45; cursor: default; }
button.ghost { background: #fff; color: #1d3557; }
.status { font-size: .95rem; font-weight: 600; min-height: 1.3em; margin-bottom: .3rem; }
.status.ok { color: #0a7d33; }
.status.bad { color: #c92f3c; }
.status.info { color: #1d3557; }
.formula { font-size: .83rem; color: #555; min-height: 1.2em; font-family: ui-monospace, monospace; }

/* Cholesky decomposition demo — build L column by column from A = LL^T */
const i18nFmt = (s, v) => s.replace(/%(\w+)%/g, (_, k) => String(v[k] ?? ''));

/* Default symmetric positive-definite 3x3 matrix A */
const DEFAULT_A = [[4, 2, -2], [2, 5, 3], [-2, 3, 14]];
let A = DEFAULT_A.map(r => r.slice());
let L = [[0,0,0],[0,0,0],[0,0,0]];
const N = 3;

/* Current position in the column-by-column sweep */
let stepJ = 0; // current column being filled
let stepI = 0; // current row within that column
let done = false;

/* Render a 3x3 matrix into the given table element */
function renderMatrix(tableId, data, highlights, emptyBelow) {
  const tbl = document.getElementById(tableId);
  tbl.innerHTML = '';
  for (let i = 0; i < N; i++) {
    const tr = document.createElement('tr');
    for (let j = 0; j < N; j++) {
      const td = document.createElement('td');
      const key = i + ',' + j;
      if (emptyBelow && j > i) {
        /* upper triangle of L is always zero */
        td.textContent = '0';
        td.className = 'zero';
      } else if (highlights && highlights.active && highlights.active === key) {
        td.textContent = fmt4(data[i][j]);
        td.className = 'active';
      } else if (highlights && highlights.done && highlights.done.has(key)) {
        td.textContent = fmt4(data[i][j]);
        td.className = 'done';
      } else if (emptyBelow && i > 0 && data[i][j] === 0 && !isDone(i, j)) {
        td.textContent = '?';
        td.className = 'empty';
      } else {
        td.textContent = fmt4(data[i][j]);
        td.className = data[i][j] === 0 ? 'zero' : '';
      }
      tr.appendChild(td);
    }
    tbl.appendChild(tr);
  }
}

function isDone(i, j) {
  /* Returns true if entry (i,j) has already been computed */
  if (j < stepJ) return true;
  if (j === stepJ && i < stepI) return true;
  return false;
}

function fmt4(v) {
  /* Format a number to at most 3 decimal places, stripping trailing zeros */
  if (v === 0) return '0';
  const s = v.toFixed(3).replace(/\.?0+$/, '');
  return s;
}

function getDoneSet() {
  /* Build the set of already-computed L entries */
  const s = new Set();
  for (let j = 0; j < stepJ; j++)
    for (let i = j; i < N; i++) s.add(i + ',' + j);
  for (let i = stepJ; i < stepI; i++) s.add(i + ',' + stepJ);
  return s;
}

function renderAll() {
  const doneSet = getDoneSet();
  const activeKey = (stepJ < N && !done) ? stepI + ',' + stepJ : null;
  renderMatrix('matA', A, null, false);
  renderMatrix('matL', L, { active: activeKey, done: doneSet }, true);
}

/* Compute the next entry of L using the Cholesky recurrence */
function choleskyStep() {
  if (done) return;
  if (stepJ >= N) { done = true; return; }

  const j = stepJ;
  const i = stepI;

  if (i === j) {
    /* Diagonal entry: l_jj = sqrt(a_jj - sum_{k<j} l_jk^2) */
    let sum = 0;
    for (let k = 0; k < j; k++) sum += L[j][k] * L[j][k];
    const val = A[j][j] - sum;
    if (val <= 0) {
      setStatus('Matrix is not positive-definite — square root of a non-positive number.', 'bad');
      done = true;
      return;
    }
    L[j][j] = Math.sqrt(val);
    setFormula(i18nFmt('l[%i%][%i%] = sqrt(%val% - %sum%) = %res%', { i: j, val: A[j][j].toFixed(3), sum: sum.toFixed(3), res: L[j][j].toFixed(4) }));
  } else {
    /* Off-diagonal entry: l_ij = (a_ij - sum_{k<j} l_ik * l_jk) / l_jj */
    let sum = 0;
    for (let k = 0; k < j; k++) sum += L[i][k] * L[j][k];
    L[i][j] = (A[i][j] - sum) / L[j][j];
    setFormula(i18nFmt('l[%i%][%j%] = (%aij% - %sum%) / %ljj% = %res%', { i: i, j: j, aij: A[i][j].toFixed(3), sum: sum.toFixed(3), ljj: L[j][j].toFixed(4), res: L[i][j].toFixed(4) }));
  }

  /* Move to the next (i,j) in column-major order within the lower triangle */
  stepI++;
  if (stepI >= N) {
    stepJ++;
    stepI = stepJ;
  }
  if (stepJ >= N) {
    done = true;
    setStatus('Done! L is complete. Verify: L × Lᵀ = A.', 'ok');
    document.getElementById('btnStep').disabled = true;
    document.getElementById('btnAll').disabled = true;
  } else {
    setStatus(i18nFmt('Computing L[%i%][%j%]…', { j: stepJ, i: stepI }), 'info');
  }
}

function setStatus(t, cls) {
  const s = document.getElementById('status');
  s.textContent = t; s.className = 'status ' + cls;
}
function setFormula(t) {
  document.getElementById('formula').textContent = t;
}

/* Reset L and step counters to restart the factorization */
function reset() {
  L = [[0,0,0],[0,0,0],[0,0,0]];
  stepJ = 0; stepI = 0; done = false;
  document.getElementById('btnStep').disabled = false;
  document.getElementById('btnAll').disabled = false;
  setStatus('Ready — click Next entry to begin.', 'info');
  setFormula('');
  renderAll();
}

/* Make the cells of A editable so the user can try different matrices */
function makeAEditable() {
  const tbl = document.getElementById('matA');
  for (let i = 0; i < N; i++) {
    for (let j = 0; j < N; j++) {
      const td = tbl.rows[i].cells[j];
      td.contentEditable = 'true';
      td.className += ' editable';
      td.addEventListener('blur', () => {
        const v = parseFloat(td.textContent);
        if (!isNaN(v)) {
          A[i][j] = v;
          A[j][i] = v; // keep A symmetric when the user edits an entry
          reset();
        } else {
          td.textContent = fmt4(A[i][j]);
        }
      });
    }
  }
}

document.getElementById('btnStep').addEventListener('click', () => { choleskyStep(); renderAll(); });
document.getElementById('btnAll').addEventListener('click', () => {
  while (!done) choleskyStep();
  renderAll();
});
document.getElementById('btnReset').addEventListener('click', () => { A = DEFAULT_A.map(r => r.slice()); reset(); makeAEditable(); });

reset();
makeAEditable();
Result
Console

Notice the pattern: the diagonal entry ljjl_{jj} is a square root (it must be real because AA is positive-definite), and each off-diagonal entry lijl_{ij} is just a dot product divided by ljjl_{jj}. This column-by-column sweep visits each entry of LL exactly once — roughly half the work of a full LU decomposition on the same matrix.

The Real Complexity

How fast is Cholesky, and why does symmetry help so much?

  • Operation count. A general LU decomposition of an n×nn \times n matrix costs about 23n3\frac{2}{3}n^3 floating-point operations. Cholesky costs roughly 13n3\frac{1}{3}n^3 — exactly half — because it only touches the lower triangle.
  • Numerical stability. For symmetric positive-definite matrices, Cholesky is unconditionally stable: no pivoting is required. The diagonal entries ljj=ajj−∑k=1j−1ljk2l_{jj} = \sqrt{a_{jj} - \sum_{k=1}^{j-1} l_{jk}^2} are always real and positive, so the algorithm never divides by zero or amplifies rounding errors the way partial pivoting must guard against in general LU.
  • Existence and uniqueness. A matrix has a Cholesky factorization if and only if it is symmetric positive-definite — and when it exists, the factorization with positive diagonal is unique. This equivalence makes "does Cholesky succeed?" a practical test for positive-definiteness.
  • Blocked and parallel variants. Modern LAPACK implementations use blocked Cholesky (DPOTRF), processing panels of columns at once to exploit cache locality and BLAS-3 matrix-matrix operations. On multi-core and GPU hardware, the triangular solves that follow (DPOTRS) are also highly parallel.

Cholesky is a solved, polished algorithm — not an open problem. Its complexity class is the same as matrix multiplication in general, and the O(n3)O(n^3) bound has been the practical standard since the 1960s. Current research focuses on sparse Cholesky (exploiting sparsity patterns in finite-element meshes) and randomized Cholesky for approximate factorizations of very large matrices.

Where It Matters

Symmetric positive-definite matrices appear wherever you measure covariance, minimize a quadratic objective, or discretize a self-adjoint differential operator — and Cholesky is almost always the right tool:

  • Statistics and machine learning. Covariance matrices are symmetric positive-semidefinite by construction. Cholesky gives you the "square root" of a covariance matrix: if Σ=LLT\Sigma = L L^T and z∼N(0,I)z \sim \mathcal{N}(0, I), then Lz∼N(0,Σ)Lz \sim \mathcal{N}(0, \Sigma). Gaussian process regression relies on this at every prediction step.
  • Kalman filtering. The "square-root Kalman filter" propagates LL (the Cholesky factor of the covariance matrix) rather than the matrix itself, giving better numerical conditioning for navigation, tracking, and control systems.
  • Convex optimization. Interior-point methods for semidefinite programming and quadratic programs solve systems Ax=bAx = b with symmetric positive-definite AA at every iteration — Cholesky is the inner loop.
  • Finite-element analysis. Discretizing an elliptic PDE produces a large, sparse, symmetric positive-definite stiffness matrix. Sparse Cholesky (with fill-reducing reordering like AMD or nested dissection) is the standard direct solver.
  • Monte Carlo simulation. Generating correlated random variables from a multivariate distribution requires factoring the correlation matrix — a textbook Cholesky application.

In short: wherever symmetry and positive-definiteness arise together — which is most of applied mathematics — Cholesky delivers the fastest, numerically cleanest linear solve available.

Conclusion

Cholesky decomposition is a beautiful case of structure meeting efficiency: when you know a matrix is symmetric and positive-definite, you don't have to pretend it isn't. Exploit that structure and the work literally halves, pivoting disappears, and numerical stability comes for free.

André-Louis Cholesky discovered this while fitting survey measurements together on a Cretan hillside. A century later, his algorithm lives inside every statistical package, every Kalman filter, every interior-point optimizer — quietly doing its job in 13n3\frac{1}{3}n^3 flops with no fuss.

The lesson generalizes: the best algorithms are not always the most general ones. A method that knows something extra about its input — symmetry, sparsity, positivity — can be twice as fast, twice as stable, and far more elegant. Understanding when those assumptions hold is as important as understanding the algorithm itself.

Share this article

Pick a channel — or use your device's native share sheet.

Comments

Loading comments...

https://www.kipuhub.com/en/article/cholesky-decomposition/Content licensed under CC BY-NC 4.0.