Skip to main content

Subsection 6.4.1 Numerical stability of triangular solve

We now use the error results for the dot product to derive a backward error result for solving \(L x = y \text{,}\) where \(L \) is an \(n \times n \) lower triangular matrix, via the algorithm in Figure 6.4.1.1, a variation on the algorithm in Figure 5.3.5.1 that stores the result in vector \(x \) and does not assume that \(L \) is unit lower triangular.

\begin{equation*} \newcommand{\routinename}{\mbox{Solve } L x = y } \newcommand{\guard}{ n( L_{TL} ) \lt n( L ) } \newcommand{\partitionings}{ L \rightarrow \FlaTwoByTwo{L_{TL}}{L_{TR}} {L_{BL}}{L_{BR}}, x \rightarrow \FlaTwoByOne{x_T}{x_B}, y \rightarrow \FlaTwoByOne{y_T}{y_B} } \newcommand{\partitionsizes}{ L_{TL} {\rm ~is~} 0 \times 0 \mbox{ and } y_T, x_T \mbox{ have } 0 \mbox{ elements} } \newcommand{\repartitionings}{ \FlaTwoByTwo{L_{TL}}{L_{TR}} {L_{BL}}{L_{BR}} \rightarrow \FlaThreeByThreeBR{L_{00}}{l_{01}}{L_{02}} {l_{10}^T}{\lambda_{11}}{l_{12}^T} {L_{20}}{l_{21}}{L_{22}}, \FlaTwoByOne{x_T}{x_B} \rightarrow \FlaThreeByOneB{x_0}{\chi_1}{x_2}, \FlaTwoByOne{y_T}{y_B} \rightarrow \FlaThreeByOneB{y_0}{\psi_1}{y_2} } \newcommand{\moveboundaries}{ \FlaTwoByTwo{L_{TL}}{L_{TR}} {L_{BL}}{L_{BR}} \leftarrow \FlaThreeByThreeTL{L_{00}}{l_{01}}{L_{02}} {l_{10}^T}{\lambda_{11}}{l_{12}^T} {L_{20}}{l_{21}}{L_{22}}, \FlaTwoByOne{x_T}{x_B} \leftarrow \FlaThreeByOneT{x_0}{\chi_1}{x_2}, \FlaTwoByOne{y_T}{y_B} \leftarrow \FlaThreeByOneT{y_0}{\psi_1}{y_2} } \newcommand{\update}{ \begin{array}{l} \chi_1 := ( \psi_1 - l_{10}^T x_0 ) / \lambda_{11} \end{array} } \FlaAlgorithm \end{equation*}
Figure 6.4.1.1. Dot product based lower triangular solve algorithm.

To establish the backward error result for this algorithm, we need to understand the error incurred in the key computation

\begin{equation*} \chi_1 := ( \psi_1 - l_{10}^T x_0 ) / \lambda_{11}. \end{equation*}

The following theorem gives the required (forward error) result, abstracted away from the specifics of how it occurs in the lower triangular solve algorithm.

Homework 6.4.1.1.

Prove Lemma 6.4.1.2

Hint

Use the Alternative Computations Model (Subsubsection 6.2.3.3) appropriately.

Solution

We know that

  • From Corollary 6.3.3.2 R-1B: if \(\beta = x^T y \) then \(\check \beta = (x + \delta\!x)^T y \) where \(\vert \delta\!x \vert \leq \gamma_n \vert x \vert \text{.}\)

  • From the ACM (Subsubsection 6.2.3.3): If \(\nu = ( \alpha - \beta ) / \lambda \) then

    \begin{equation*} \check \nu = \frac{ \alpha - \beta }{\lambda} \frac{1}{( 1 + \epsilon_- ) (1 + \epsilon_/)}, \end{equation*}

    where \(\vert \epsilon_-\vert \leq \meps \) and \(\vert \epsilon_/\vert \leq \meps \text{.}\)

Hence

\begin{equation*} \check \nu = \frac{\alpha - ( x + \delta\!x)^T y }{\lambda } \frac{1}{( 1 + \epsilon_{-} ) (1 + \epsilon_{/})}, \end{equation*}

or, equivalently,

\begin{equation*} \lambda ( 1 +\epsilon_- ) (1 + \epsilon_/) \check \nu = \alpha - ( x + \delta\! x )^T y , \end{equation*}

or,

\begin{equation*} \lambda ( 1 + \theta_2 ) \check \nu = \alpha - ( x + \delta\! x )^T y , \end{equation*}

where \(\vert \theta_2 \vert \leq \gamma_2 \text{,}\) which can also be written as

\begin{equation*} ( \lambda + \delta\!\lambda ) \check \nu = \alpha - ( x + \delta\!x)^T y , \end{equation*}

where \(\delta\!\lambda = \theta_2 \lambda \) and hence \(\vert \delta\!\lambda \vert \leq \gamma_2 \| \lambda \vert \text{.}\)

The errror result for the algorithm in Figure 6.4.1.1 is given by

The reasoning behind the result is that one expects the maximal error to be incurred during the final iteration when computing \(\chi_1 := ( \psi_1 - l_{10}^T x_0 ) / \lambda_{11} \text{.}\) This fits Lemma 6.4.1.2, except that this assignment involves a dot product with vectors of length \(n-1 \) rather than of length \(n \text{.}\)

You now prove Theorem 6.4.1.3 by first proving the special cases where \(n = 1 \) and \(n = 2 \text{,}\) and then the general case.

Homework 6.4.1.2.

Prove Theorem 6.4.1.3 for the case where \(n = 1 \text{.}\)

Solution

Case 1: \(n =1 \text{.}\)

The system looks like \(\lambda_{11} \chi_{1} = \psi_{1} \) so that

\begin{equation*} \chi_1 = \psi_1 / \lambda_{11} \end{equation*}

and

\begin{equation*} \check \chi_1 = \psi_1 / \lambda_{11} \frac{1}{1 + \epsilon_{/}} \end{equation*}

Rearranging gives us

\begin{equation*} \lambda_{11} \check \chi_1 (1 + \epsilon_{/}) = \psi_1 \end{equation*}

or

\begin{equation*} ( \lambda_{11} + \delta\!\lambda_{11} ) \check \chi_1 = \psi_1 \end{equation*}

where \(\delta\!\lambda_{11} = \epsilon_{/} \lambda_{11} \) and hence

\begin{equation*} \begin{array}{rcl} \vert \delta\!\lambda_{11} \vert \amp = \amp \vert \epsilon_{/} \vert \vert \lambda_{11} \vert \\ \amp \leq \amp \gamma_1 \vert \lambda_{11} \vert \\ \amp \leq \amp \gamma_2 \vert \lambda_{11} \vert \\ \amp \leq \amp \max( \gamma_2, \gamma_{n-1} ) \vert \lambda_{11} \vert . \end{array} \end{equation*}
Homework 6.4.1.3.

Prove Theorem 6.4.1.3 for the case where \(n = 2 \text{.}\)

Solution

Case 2: \(n =2 \text{.}\)

The system now looks like

\begin{equation*} \left( \begin{array}{c | c} \lambda_{00} \amp 0 \\ \hline \lambda_{10} \amp \lambda_{11} \end{array} \right) \left( \begin{array}{c} \chi_0 \\ \chi_1 \end{array} \right) = \left( \begin{array}{c} \psi_0 \\ \psi_1 \end{array} \right). \end{equation*}

From the proof of Case 1 we know that

\begin{equation} (\lambda_{00} + \delta\!\lambda_{00}) \check \chi_0 = \psi_{0}, \mbox{ where } \vert \delta\!\lambda_{00} \vert \leq \gamma_1 \vert \lambda_{00} \vert .\label{trsv-eq-1}\tag{6.4.1} \end{equation}

Since \(\chi_1 = ( \psi_1 - \lambda_{10} \check \chi_0 )/ \lambda_{11} \text{,}\) Lemma 6.4.1.2 tells us that

\begin{equation} ( \lambda_{10} + \delta\!\lambda_{10} ) \check \chi_0 + ( \lambda_{11} +\delta\!\lambda_{11} ) \check \chi_1 = \psi_1,\label{trsv-eq-2}\tag{6.4.2} \end{equation}

where

\begin{equation*} \vert \delta\!\lambda_{10} \vert \leq \gamma_1\vert \lambda_{10} \vert \mbox{ and } \vert \delta\!\lambda_{11} \vert \leq \gamma_2 \vert \lambda_{11} \vert . \end{equation*}

(6.4.1) and (6.4.2) can be combined into

\begin{equation*} \left( \begin{array}{c | c} \lambda_{00} + \delta\!\lambda_{00} \amp 0 \\ \hline \lambda_{10} + \delta\!\lambda_{10} \amp \lambda_{11} + \delta\!\lambda_{11} \end{array} \right) \left( \begin{array}{c} \check \chi_0 \\ \check \chi_1 \end{array} \right) = \left( \begin{array}{c} \psi_0 \\ \psi_1 \end{array} \right), \end{equation*}

where

\begin{equation*} \left( \begin{array}{c | c} \vert \delta\!\lambda_{00} \vert \amp 0 \\ \hline \vert \delta\!\lambda_{10} \vert \amp \vert \delta\!\lambda_{11} \vert \end{array} \right) \leq \left( \begin{array}{c | c} \gamma_1 \vert \lambda_{00} \vert\amp 0 \\ \hline \gamma_1 \vert \lambda_{10} \vert \amp \gamma_2 \vert \lambda_{11} \vert \end{array} \right). \end{equation*}

Since \(\gamma_1 \leq \gamma_2 \)

\begin{equation*} \left\vert \left( \begin{array}{c | c} \delta\!\lambda_{00} \amp 0 \\ \hline \delta\!\lambda_{10} \amp \delta\!\lambda_{11} \end{array} \right) \right\vert \leq \max( \gamma_2, \gamma_{n-1} ) \left\vert \left( \begin{array}{c | c} \lambda_{00} \amp 0 \\ \hline \lambda_{10} \amp \lambda_{11} \end{array} \right) \right\vert. \end{equation*}
Homework 6.4.1.4.

Prove Theorem 6.4.1.3 for \(n \geq 1 \text{.}\)

Solution

We will utilize a proof by induction.

  • Case 1: \(n =1 \text{.}\)

    See Homework 6.4.1.2.

  • Case 2: \(n =2 \text{.}\)

    See Homework 6.4.1.3.

  • Case 3: \(n \gt 2 \text{.}\)

    The system now looks like

    \begin{equation} \left( \begin{array}{c | c} L_{00} \amp 0 \\ \hline l_{10}^T \amp \lambda_{11} \end{array} \right) \left( \begin{array}{c} x_0 \\ \chi_1 \end{array} \right) = \left( \begin{array}{c} y_0 \\ \psi_1 \end{array} \right) ,\label{trsv-eq-1c}\tag{6.4.3} \end{equation}

    where \(L_{00} \in \R^{(n-1) \times (n-1)} \text{,}\) and the inductive hypothesis states that

    \begin{equation*} (L_{00} + \Delta\!L_{00}) \check x_0 = y_0 \mbox{ where } \vert \Delta\!L_{00} \vert \leq \max( \gamma_2, \gamma_{n-2} ) \vert L_{00} \vert. \end{equation*}

    Since \(\chi_1 = ( \psi_1 - l_{10}^T \check x_0 )/ \lambda_{11} \text{,}\) Lemma 6.4.1.2 tells us that

    \begin{equation} ( l_{10} + \delta\!l_{10})^T\check x_0 + (\lambda_{11} + \delta\!\lambda_{11}) \check \chi_1 = \psi_1,\label{trsv-eq-2c}\tag{6.4.4} \end{equation}

    where \(\vert \delta\!l_{10} \vert \leq \gamma_{n-1} \vert l_{10} \vert \) and \(\vert \delta\!\lambda_{11} \vert \leq \gamma_2 \vert \lambda_{11} \vert \text{.}\)

    (6.4.3) and (6.4.4) can be combined into

    \begin{equation*} \left( \begin{array}{c | c} L_{00} + \delta\!L_{00} \amp 0 \\ \hline (l_{10} + \delta\!l_{10})^T \amp \lambda_{11} + \delta\!\lambda_{11} \end{array} \right) \left( \begin{array}{c} \check x_0 \\ \check \chi_1 \end{array} \right) = \left( \begin{array}{c} y_0 \\ \psi_1 \end{array} \right), \end{equation*}

    where

    \begin{equation*} \left( \begin{array}{c | c} \vert \delta\!L_{00} \vert \amp 0 \\ \hline \vert \delta\!l_{10}^T \vert \amp \vert \delta\!\lambda_{11} \vert \end{array} \right) \leq \left( \begin{array}{c | c} \max( \gamma_2, \gamma_{n-2} ) \vert L_{00} \vert\amp 0 \\ \hline \gamma_{n-1} \vert l_{10}^T \vert \amp \gamma_2 \vert \lambda_{11} \vert \end{array} \right) \end{equation*}

    and hence

    \begin{equation*} \left\vert \left( \begin{array}{c | c} \delta\!L_{00} \amp 0 \\ \hline \delta\!l_{10}^T \amp \delta\!\lambda_{11} \end{array} \right) \right\vert| \leq \max( \gamma_2, \gamma_{n-1} ) \left\vert \left( \begin{array}{c | c} L_{00} \amp 0 \\ \hline l_{10}^T \amp \lambda_{11} \end{array} \right) \right\vert. \end{equation*}
  • By the Principle of Mathematical Induction, the result holds for all \(n \geq 1 \text{.}\)

A careful examination of the solution to Homework 6.4.1.2, together with the fact that \(\gamma_{n-1} \leq \gamma_n \) allows us to state a slightly looser, but cleaner, result of Theorem 6.4.1.3: