a really simple iterative refinement example
- 2 minutes read - 340 wordsI wanted to have a quick example of what iterative refinement is and looks like.
Suppose you wanted to solve the (very simple!) problem of solving 5.1x=16, without dividing by decimal numbers. (Yes, this is a bit contrived–I’m not sure why you would be able to divide by an integer to get a decimal, but not divide by a decimal to get a decimal. Bear with me.)
The standard way would be, of course, to simply evaluate x = (5.1)^{-1}16=3.137254. But this is not available–you can’t find an inverse of multiplying by 5.1!
So instead, you start with some point x_0=0. You calculate the residual, r_0 = 16-5.1cdot x_0=16. (Multiplying decimals is allowed. Just not inverses.) Now you solve the simpler problem 5d_0=r_0 for d_0, or 5d_0=16. Because now the coefficient is integral, you get d_0=(5)^{-1}16=3.2. Then we generate a new update for x via x_1=x_0+d_0 = 3.2.
This is already somewhat close!
Can we make it better? Well, sure. Just repeat the process:
r_1 = 16-5.1*3.2=-.32 d_1 = (5)^{-1}r_1 = (5)^{-1}(-.32) = -.064 x_2 = x_1 + d_1 = 3.2 - 0.064 = 3.136 Repeating the process a few more times:
| iteration | iterate | residual | correction |
|---|---|---|---|
| 0 | 0.000000000e+00 | 1.600000000e+01 | 3.200000000e+00 |
| 1 | 3.200000000e+00 | -3.200000000e-01 | -6.400000000e-02 |
| 2 | 3.136000000e+00 | 6.400000000e-03 | 1.280000000e-03 |
| 3 | 3.137280000e+00 | -1.280000000e-04 | -2.560000000e-05 |
| 4 | 3.137254400e+00 | 2.559999997e-06 | 5.119999994e-07 |
| 5 | 3.137254912e+00 | -5.119999713e-08 | -1.023999943e-08 |
| The algorithm for this looks like this: Suppose: |
you want to solve alpha x = beta, you can’t calculate alpha^{-1}, you can calculate hat{alpha}^{-1}, hat{alpha}approx alpha. Start with x_0 given and let i=0.
Calculate the residual r_i=beta - alpha cdot x_i. Calculate d_i=hat{alpha}^{-1} r_i. Finally, update the iterate: x_{i+1} = x_i + d_i. Now then, here’s the punchline. What if alpha is a matrix? Then it turns out that alpha^{-1} isn’t impossible to compute, just expensive. And also that hat{alpha}^{-1} is simply much faster to compute. Everything else we’ve said still applies–we can still use iterative refinement until the problem is solved.
There are some bounds for when we would expect this to work. But I’ll address that another day.