Reusing LU factors for multiple right-hand sides (source code)

= Reusing LU factors for multiple right-hand sides
{c}
{title2=$AX=B,\quad LY=PB,\quad UX=Y$}

When a fixed <matrix> $A$ is used with several forcing or measurement <vectors>, compute its <Pivoted LU decomposition> once. Collect $r$ right-hand sides into the columns of $B=(b^{(1)}\ \cdots\ b^{(r)})$, and the corresponding solutions into $X$. Then
$$
AX=B,\qquad LY=PB,\qquad UX=Y.
$$
Each column undergoes the same row permutation, <forward substitution> and <back substitution>. Neither $L$, $U$ nor $P$ depends on $B$. This is the useful separation between preparing a <linear system> and applying its solution procedure repeatedly. LAPACK implements it by taking previously computed factors, permuting all right-hand sides, and performing the two triangular solves.<DGETRS (LAPACK 3.3.1) by LAPACK project>[LAPACK DGETRS, version 3.3.1, April 2011, Purpose, NRHS and B arguments, and nontransposed branch: DLASWP followed by lower-unit and upper-nonunit DTRSM calls]{ref}

For the unpivoted factors in <Computed three-by-three LU decomposition>, take
$$
B=\begin{pmatrix}7&0&4\\23&2&6\\26&4&-8\end{pmatrix},\qquad Y=\begin{pmatrix}7&0&4\\9&2&-2\\6&-2&2\end{pmatrix},\qquad X=\begin{pmatrix}1&0&2\\2&1&-1\\3&-1&1\end{pmatrix}.
$$
Direct <matrix multiplication> checks all three identities $LY=B$, $UX=Y$ and $AX=B$. The first column reproduces the previous solve; the other two use exactly the same factors. New right-hand sides can also arrive one at a time. If $A$ changes, these stored factors generally need to be recomputed.

\Image[/lu-pivot-repeated-solves.png]
{title=One row permutation fixes a zero <pivot element>; each new right-hand side reuses the same <LU decomposition> and two triangular solves}