Numerical Analysis/Iterative Refinement

What is Iterative Refinement?

A detailed discussion of iterative refinement can be found on the Wikipedia page.

To give a brief description, it is a technique used to improve the approximate solution $\hat{x}$ to linear system $Ax=b.$ This technique is generally only used on systems that are thought or determined to be ill-conditioned.

The process involves three primary steps:

Iterative Refinement Process
Once an approximation to the solution, ${\hat{x}}_1$, has been made, with either rounded-Gaussian elimination or otherwise, complete the following steps: Compute residual: $r_i=b-A{\hat{x}}_i$ Solve $A{y}_i=r_i$ Update ${\hat{x}}_{i+1}={\hat{x}}_i+y_i$ Continue to update step one with the newly calculated ${\hat{x}}_{i+1}$ until a desired tolerance has been reached.

Approximating the Condition Number

In theory, the condition number of a matrix depends on the norms of the matrix and its inverse. In practice, however, calculating the inverse is subject to round-off error and the accuracy of the calculations.

If these calculations involve arithmetic with t digits of accuracy, then the approximate condition number for the matrix - call it A - is the norm of A the norm of the approximation of the inverse of A.

Assuming that the approximate solution to the linear system is being determined with t-digit arithmetic and Gaussian elimination, we can show

$\Vert \begin{array}{c}r\end{array}\Vert \approx 10^{-t}\Vert \begin{array}{c}A\end{array}\Vert \Vert \begin{array}{c}\hat{x}\end{array}\Vert ,$

where r is the residual vector for the approximation $\hat{x}.$

This approximation for the condition number can be obtained without having to invert matrix A.

When doing iterative refinement problems, we will have, or will calculate the values for $A,y,$ and $\hat{x}$. The approximate solution, $\hat{y}$ satisfies

$\hat{y}\approx A^{-1}r=A^{-1}(b-A\hat{x})=A^{-1}b-A^{-1}A\hat{x}=x-\hat{x},$

so $\hat{y}$ is an estimate of the error for approximate solution $\hat{x}.$ Observe that

$\begin{array}{cc}\Vert \begin{array}{c}\hat{y}\end{array}\Vert \approx \Vert \begin{array}{c}x-\hat{x}\end{array}\Vert & =\Vert \begin{array}{c}{A}^{-1}r\end{array}\Vert \le \Vert \begin{array}{c}{A}^{-1}\end{array}\Vert \Vert \begin{array}{c}r\end{array}\Vert \approx \Vert \begin{array}{c}{A}^{-1}\end{array}\Vert (10^{-t}\Vert \begin{array}{c}A\end{array}\Vert \Vert \begin{array}{c}\hat{x}\end{array}\Vert )\\ & =10^{-t}\Vert \begin{array}{c}\hat{x}\end{array}\Vert K(A).\end{array}$

This approximates the condition number associated with solving $Ax=b$ and can be re-written and expressed as seen here:

$K(A)\approx \frac{\Vert \begin{array}{c}\hat{y}\end{array}\Vert }{\Vert \begin{array}{c}\hat{x}\end{array}\Vert }10^t.$

Example

Consider the following ill-conditioned linear system

$\begin{array}{cc}3.9x_1+1.6x_2& =5.5\\ 6.8x_1+2.9x_2& =9.7\end{array}$

Part 1

Solve using Gaussian elimination and iterative refinement (and using two-digit rounding arithmetic). Continue using iterative refinement until ${\hat{x}}_4$ is found. Compare with the true solution.

Solution:

First, use a calculator or computer to find the true solution $x$ of this system. We will discover that $x=\left[\begin{array}{c}1\\ 1\end{array}\right].$

Set up the augmented matrix and solve using Guassian elimination, making sure to use only two-digit rounding arithmetic.

$\begin{array}{ccc}\left[\begin{array}{ccc}3.9& 1.6& 5.5\\ 6.8& 2.9& 9.7\end{array}\right]& =\begin{array}{c}R1/3.9\\ R2/6.8\end{array}& \to \left[\begin{array}{ccc}1& 0.41& 1.4\\ 1& 0.43& 1.4\end{array}\right]\\ \\ & =\begin{array}{c}R2-R1\end{array}& \to \left[\begin{array}{ccc}1& 0.41& 1.4\\ 0& 0.02& 0\end{array}\right]\\ \\ & =\begin{array}{c}R2/0.02\end{array}& \to \left[\begin{array}{ccc}1& 0.41& 1.4\\ 0& 1& 0\end{array}\right]\end{array}$

Solving with back-substitution we find the first approximation to the system:

${\hat{x}}_1=\left[\begin{array}{c}{x}_1\\ x_2\end{array}\right]=\left[\begin{array}{c}1.4\\ 0\end{array}\right]$

Because we mentioned that this system is ill-conditioned, we can use the iterative refinement process to find a better approximation. Recall the steps mentioned above:

1. Compute residual: $r_i=b-A{\hat{x}}_i$
   
2. Solve $A{y}_i=r_i$
   
3. Update ${\hat{x}}_{i+1}={\hat{x}}_i+y_i$
   
   

STEP 1: Compute the residual $r_1.$

$\begin{array}{cc}{r}_1=b-A{\hat{x}}_1& =\left[\begin{array}{c}5.5\\ 9.7\end{array}\right]-\left[\begin{array}{cc}3.9& 1.6\\ 6.8& 2.9\end{array}\right]\left[\begin{array}{c}1.4\\ 0\end{array}\right]\\ \\ & =\left[\begin{array}{c}5.5\\ 9.7\end{array}\right]-\left[\begin{array}{c}5.5\\ 9.5\end{array}\right]\\ \\ & =\left[\begin{array}{c}0\\ 0.2\end{array}\right]=r_1\end{array}$

STEP 2: Now solve $A{y}_1=r_1$ using rounded Gaussian elimination to find the $y$ vector.

$\begin{array}{cc}A{y}_1=r_1\Rightarrow \left[\begin{array}{cc}3.9& 1.6\\ 6.8& 2.9\end{array}\right]\left[\begin{array}{c}{y}_1\\ y_2\end{array}\right]=\left[\begin{array}{c}0\\ 0.2\end{array}\right]& \Rightarrow \left[\begin{array}{ccc}3.9& 1.6& 0\\ 6.8& 2.9& 0.2\end{array}\right]=\left[\begin{array}{ccc}1& 0.41& 0\\ 1& 0.43& 0.029\end{array}\right]\\ \\ & =\left[\begin{array}{ccc}1& 0.41& 0\\ 0& 0.02& 0.029\end{array}\right]=\left[\begin{array}{ccc}1& 0.41& 0\\ 0& 1& 1.5\end{array}\right]\\ \\ & =\left[\begin{array}{ccc}1& 0& -0.62\\ 0& 1& 1.5\end{array}\right]\\ \\ \end{array}$

$\Rightarrow y_1=\left[\begin{array}{c}-0.62\\ 1.5\end{array}\right]$

STEP 3: Update $\hat{x}.$

$\begin{array}{cc}{\hat{x}}_2& ={\hat{x}}_1+y_1\\ \\ & =\left[\begin{array}{c}1.4\\ 0\end{array}\right]+\left[\begin{array}{c}-0.62\\ 1.5\end{array}\right]\\ \\ & =\left[\begin{array}{c}0.78\\ 1.5\end{array}\right]\end{array}$

Note that ${\hat{x}}_2$ is a better representation of $x$, than ${\hat{x}}_1,$ however it is still not a very accurate representation. Continue with a couple more iterations.

Repeat, as needed.

In the second iteration of Iterative Refinement, we find the following:

$r_2=\left[\begin{array}{c}0.1\\ 0\end{array}\right],y_2=\left[\begin{array}{c}0.56\\ -1.3\end{array}\right]$

Thus,

$\begin{array}{cc}{\hat{x}}_3& ={\hat{x}}_2+y_2=\left[\begin{array}{c}0.78\\ 1.5\end{array}\right]+\left[\begin{array}{c}0.56\\ -1.3\end{array}\right]=\left[\begin{array}{c}1.3\\ 0.2\end{array}\right]\end{array}$

In the third iteration we find

$r_3=\left[\begin{array}{c}0.1\\ -0.3\end{array}\right],y_3=\left[\begin{array}{c}-0.34\\ 0.9\end{array}\right],$

which we can use to find our final approximation, ${\hat{x}}_4$ as

$\begin{array}{cc}{\hat{x}}_4& ={\hat{x}}_3+y_3=\left[\begin{array}{c}1.3\\ 0.2\end{array}\right]+\left[\begin{array}{c}-0.34\\ 0.9\end{array}\right]=\left[\begin{array}{c}0.96\\ 1.1\end{array}\right]\end{array}$

Recall that the true solution of the linear system is $x=\left[\begin{array}{c}1\\ 1\end{array}\right].$ After completing three iterations of iterative refinement we have approximate the solution ${\hat{x}}_4=\left[\begin{array}{c}0.96\\ 1.1\end{array}\right],$ which is significantly closer than our original approximation ${\hat{x}}_1=\left[\begin{array}{c}1.4\\ 0\end{array}\right].$

Part 2

Estimate the condition number $K(A)$ in Part 1 using the method discussed above, and compare this with the condition number calculated using norms.

Solution:

We can find the true condition number by calculating it using norms - in the way it is described on the condition number Wikipedia . (Keep in mind, we will still use 2 significant digits in this example as well.)

$A=\left[\begin{array}{cc}3.9& 1.6\\ 6.8& 2.9\end{array}\right],$ and $A^{-1}=\left[\begin{array}{cc}6.7& -3.7\\ -16& 9.1\end{array}\right].$

We can then solve for the condition number $K(A).$

$K(A)=\Vert \begin{array}{c}A\end{array}\Vert \Vert \begin{array}{c}{A}^{-1}\end{array}\Vert =\Vert \begin{array}{c}\left[\begin{array}{cc}3.9& 1.6\\ 6.8& 2.9\end{array}\right]\end{array}\Vert \Vert \begin{array}{c}\left[\begin{array}{cc}6.7& -3.7\\ -16& 9.1\end{array}\right]\end{array}\Vert =\Vert \begin{array}{c}9.7\end{array}\Vert \Vert \begin{array}{c}25\end{array}\Vert =240.$

We can use our approximations, ${\hat{x}}_1$ and $y_1$, found in the first iteration of Example 1 to use in the method we learned above, on Approximating the Condition Number.

$\begin{array}{c}K(A)\approx \frac{\Vert \begin{array}{c}{\hat{y}}_1\end{array}\Vert }{\Vert \begin{array}{c}{\hat{x}}_1\end{array}\Vert }10^t=\frac{\Vert \begin{array}{c}\left[\begin{array}{c}-0.62\\ 1.5\end{array}\right]\end{array}\Vert }{\Vert \begin{array}{c}\left[\begin{array}{c}1.4\\ 0\end{array}\right]\end{array}\Vert }10^2=\frac{1.5}{1.4}10^2=1.07*10^2=107.\end{array}$

So we have found the true condition number to be $240$ and the approximation to be $107$. These numbers are both in the realm of $10^2$, meaning both for the true and approximate condition number, we can see there may be a loss of two digits of accuracy.

Quiz

<quiz display=simple>

Not only does iterative refinement improve ill-conditioned matrices, it is also very help in improving well-conditioned matrices. |type="()"} +TRUE -FALSE

{Which of the following matrices might benefit from using iterative refinement? |type="[]"} + $\left[\begin{array}{cc}1& 2\\ 1.0001& 2\end{array}\right]$ + $\left[\begin{array}{cc}1& 2\\ 10,001& 2\end{array}\right]$ - $\left[\begin{array}{cc}1& 2\\ 10,001& 20,002\end{array}\right]$ - $\left[\begin{array}{cc}1& 20,002\\ 10,001& 2\end{array}\right]$

{ |type="{}"} If you calculate another step of iterative refinement on Example 1, what will be the new approximation, ${\hat{x}}_5=\left[\begin{array}{c}{x}_1\\ x_2\end{array}\right]$  ? $x_1=$ { 0.96 } $x_2=$ { 1.1 }

{ |type="{}"} A linear system is given by

$\left[\begin{array}{ccc}3.3330& 15920& -10.333\\ 2.2220& 16.710& 9.6120\\ 1.5611& 5.1791& 1.6852\end{array}\right]$ $\left[\begin{array}{c}{x}_1\\ x_2\\ x_3\end{array}\right]$ = $\left[\begin{array}{c}15913\\ 28.544\\ 8.4254\end{array}\right],$ with ${\hat{x}}_1$ = $\left[\begin{array}{c}1.2001\\ 0.99991\\ 0.92538\end{array}\right],$ and $y_1$ = $\left[\begin{array}{c}-0.20008\\ 8.9987\times 10^{-5}\\ 0.074607\end{array}\right].$

Based on this info, what is the condition number using the approximation and the typical norms-multiplication method? Using Approximation $K(A)=$ { 16672 } Using Norms $K(A)=$ { 15999 }

</quiz>

References

Burden and Faires. Numerical Analysis, 3rd Edition. Template:ISBN