To add preconditioning to
\begin{equation*}
A x = b
\end{equation*}
we pick a SPD preconditioner \(M = \tilde L \tilde L^T \) and instead solve the equivalent problem
\begin{equation*}
\begin{array}[t]{c}
\underbrace{
\tilde L^{-1} A \tilde L^{-T}
}
\\
\tilde A
\end{array}
\begin{array}[t]{c}
\underbrace{
\tilde L^T x
}
\\
\tilde x
\end{array}
=
\begin{array}[t]{c}
\underbrace{
\tilde L^{-1} b .
} \\
\tilde b
\end{array}
\end{equation*}
This changes the algorithm in Figure 8.3.5.2 (right) to
\begin{equation*}
\begin{array}{l}
{\bf Given:~} A, b, M = \tilde L \tilde L^T \\
\tilde x^{(0)} := 0
\\
\tilde A = \tilde L^{-1} A \tilde L^{-T} \\
\tilde r^{(0)} := \tilde L^{-1} b \\
k := 0 \\
{\bf while~} \tilde r^{(k)} \neq 0 \\
~~~ {\bf if~} k = 0 \\
~~~ ~~~ \tilde p^{(k)} = \tilde r^{(0)} \\
~~~ {\bf else} \\
~~~ ~~~ \tilde \gamma_k := (\tilde r^{(k)\,T} \tilde r^{(k)}) /
(\tilde r^{(k-1)\,T} \tilde r^{(k-1)}) \\
~~~ ~~~ \tilde p^{(k)} := \tilde r^{(k)} + \tilde
\gamma_k \tilde p^{(k-1)}
\\
~~~ {\bf endif} \\
~~~ \tilde \alpha_k := \frac{\tilde r^{(k)\,T} \tilde
r^{(k)}}{\tilde p^{(k)\,T}
\tilde A \tilde p^{(k)}}
\\
~~~ \tilde x^{(k+1)} := \tilde x^{(k)} + \tilde
\alpha_k \tilde p^{(k)}
\\
~~~ \tilde r^{(k+1)} := \tilde r^{(k)} - \tilde
\alpha_k \tilde A \tilde p^{(k)} \\
~~~ k := k + 1 \\
{\bf endwhile}
\end{array}
\end{equation*}
Now, much like we did in the constructive solution to Homework 8.2.5.1 we now morph this into an algorithm that more directly computes \(x^{(k+1)} \text{.}\) We start by substituting
\begin{equation*}
\tilde A =\tilde L^{-1} A \tilde L^{-T},
\tilde x^{(k)} = \tilde L^T x^{(k)}, \tilde r^{(k)} = \tilde
L^{-1} r^{(k)}, \tilde p^{(k)} = \tilde L^{T} p^{(k)},
\end{equation*}
which yields
\begin{equation*}
\begin{array}{l}
{\bf Given:~} A, b, M = \tilde L \tilde L^T \\
\tilde L^T x^{(0)} := 0
\\
\tilde L^{-1} r^{(0)} := \tilde L^{-1} b \\
k := 0 \\
{\bf while~} \tilde L^{-1} r^{(k)} \neq 0 \\
~~~ {\bf if~} k = 0 \\
~~~ ~~~ \tilde L^{T} p^{(k)} = \tilde L^{-1} r^{(0)} \\
~~~ {\bf else} \\
~~~ ~~~ \tilde \gamma_k := ((\tilde L^{-1} r^{(k)})^T
\tilde L^{-1} r^{(k)}) /
(\tilde L^{-1} r^{(k-1)})^T \tilde L^{-1} r^{(k-1)}) \\
~~~ ~~~ \tilde L^{T} p^{(k)} := \tilde L^{-1} r^{(k)} + \tilde
\gamma_k \tilde L^{T} p^{(k-1)}
\\
~~~ {\bf endif} \\
~~~ \tilde \alpha_k := \frac{(\tilde L^{-1} r^{(k)})^T \tilde
L^{-1} r^{(k)}}{((\tilde L^{T} p^{(k)})^T
\tilde L^{-1} A \tilde L^{-T} \tilde L^{T} p^{(k)}}
\\
~~~ \tilde L^T x^{(k+1)} := \tilde L^T x^{(k)} + \tilde
\alpha_k \tilde L^{T} p^{(k)}
\\
~~~ \tilde L^{-1} r^{(k+1)} := \tilde L^{-1} r^{(k)} - \tilde
\alpha_k \tilde L^{-1} A \tilde L^{-T} \tilde L^{T} p^{(k)} \\
~~~ k := k + 1 \\
{\bf endwhile}
\end{array}
\end{equation*}
If we now simplify and manipulate various parts of this algorithm we get
\begin{equation*}
\begin{array}{l}
{\bf Given:~} A, b, M = \tilde L \tilde L^T \\
x^{(0)} := 0
\\
r^{(0)} := b \\
k := 0 \\
{\bf while~} r^{(k)} \neq 0 \\
~~~ {\bf if~} k = 0 \\
~~~ ~~~ p^{(k)} = M^{-1} r^{(0)} \\
~~~ {\bf else} \\
~~~ ~~~ \tilde \gamma_k := (r^{(k)\,T} M^{-1} r^{(k)}) /
(r^{(k-1)\,T} M^{-1} r^{(k-1)}) \\
~~~ ~~~ p^{(k)} := M^{-1} r^{(k)} + \tilde
\gamma_k p^{(k-1)}
\\
~~~ {\bf endif} \\
~~~ \tilde \alpha_k := \frac{r^{(k)\,T} M^{-1} r^{(k)}}{p^{(k)\,T} A p^{(k)}}
\\
~~~ x^{(k+1)} := x^{(k)} + \tilde
\alpha_k p^{(k)}
\\
~~~ r^{(k+1)} := r^{(k)} - \tilde
\alpha_k A p^{(k)} \\
~~~ k := k + 1 \\
{\bf endwhile}
\end{array}
\end{equation*}
Finally, we avoid the recomputing of \(M^{-1} r^{(k)}\) and \(A p^{(k)} \) by introducing \(z^{(k)} \) and \(q^{(k)} \text{:}\)
\begin{equation*}
\begin{array}{l}
{\bf Given:~} A, b, M = \tilde L \tilde L^T \\
x^{(0)} := 0
\\
r^{(0)} := b \\
k := 0 \\
{\bf while~} r^{(k)} \neq 0 \\
~~~ z^{(k)} := M^{-1} r^{(k)} \\
~~~ {\bf if~} k = 0 \\
~~~ ~~~ p^{(k)} = z^{(0)} \\
~~~ {\bf else} \\
~~~ ~~~ \tilde \gamma_k := (r^{(k)\,T} z^{(k)}) /
(r^{(k-1)\,T} z^{(k-1)}) \\
~~~ ~~~ p^{(k)} := z^{(k)} + \tilde
\gamma_k p^{(k-1)}
\\
~~~ {\bf endif} \\
~~~ q^{(k)} := A p^{(k)} \\
~~~ \tilde \alpha_k := \frac{r^{(k)\,T} z^{(k)}}{p^{(k)\,T} q^{(k)}}
\\
~~~ x^{(k+1)} := x^{(k)} + \tilde
\alpha_k p^{(k)}
\\
~~~ r^{(k+1)} := r^{(k)} - \tilde
\alpha_k q^{(k)} \\
~~~ k := k + 1 \\
{\bf endwhile}
\end{array}
\end{equation*}
(Obviously, there are a few other things that can be done to avoid unnecessary recomputations of \(r^{(k)\,T}
z^{(k)} \text{.}\))
