Problem

MATLAB

Consider the following program,

formatSpec = 'With %d sqrt, then %d times ^2 operations, the number %.16f becomes: %.16f \n'; % the string format for fprintf function
for n = 1:60
    r_original = 2.0;
    r = r_original;
    for i = 1:n
        r = sqrt(r);
    end
    for i = 1:n
        r = r^2;
    end
    fprintf(formatSpec,n,n,r_original,r);
end

Explain what this code does. Then run the code, and explain why do you see the behavior observed. In particular, why do you not recover the original value $2.0$ after many repetitions of the same forward and reverse task of taking the square root and squaring the result?

Python

Consider the following program,

from math import sqrt
for n in range(1, 60):
    r_org = 2.0
    r = r_org
    for i in range(n):
        r = sqrt(r)
    for i in range(n):
        r = r ** 2
    print ('With {} times sqrt and then {} times **2, the number {} becomes: {:.16f}'.format(n,n,r_org,r))

Explain what this code does. Then run the code, and explain why do you the behavior observed. In particular, why do you not recover the original value $2$ after many repetitions of the same forward and reverse task?

Solution

MATLAB

This code (roundoff.m) will yield the following output,

roundoff
With 1 sqrt, then 1 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000000000004 
With 2 sqrt, then 2 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999999999996 
With 3 sqrt, then 3 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999999999996 
With 4 sqrt, then 4 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999999999964 
With 5 sqrt, then 5 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999999999964 
With 6 sqrt, then 6 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999999999964 
With 7 sqrt, then 7 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999999999714 
With 8 sqrt, then 8 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000000000235 
With 9 sqrt, then 9 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000000000235 
With 10 sqrt, then 10 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000000000235 
With 11 sqrt, then 11 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000000000235 
With 12 sqrt, then 12 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999999991336 
With 13 sqrt, then 13 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999999973292 
With 14 sqrt, then 14 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999999973292 
With 15 sqrt, then 15 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999999973292 
With 16 sqrt, then 16 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000000117746 
With 17 sqrt, then 17 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000000408580 
With 18 sqrt, then 18 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000000408580 
With 19 sqrt, then 19 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000001573586 
With 20 sqrt, then 20 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000001573586 
With 21 sqrt, then 21 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000001573586 
With 22 sqrt, then 22 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000010885857 
With 23 sqrt, then 23 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000029511749 
With 24 sqrt, then 24 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000066771721 
With 25 sqrt, then 25 times ^2 operations, the number 2.0000000000000000 becomes: 2.0000000066771721 
With 26 sqrt, then 26 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999917775542 
With 27 sqrt, then 27 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999917775542 
With 28 sqrt, then 28 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999917775542 
With 29 sqrt, then 29 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999917775542 
With 30 sqrt, then 30 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999917775542 
With 31 sqrt, then 31 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999999917775542 
With 32 sqrt, then 32 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999990380770896 
With 33 sqrt, then 33 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999971307544144 
With 34 sqrt, then 34 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999971307544144 
With 35 sqrt, then 35 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999971307544144 
With 36 sqrt, then 36 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999971307544144 
With 37 sqrt, then 37 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999971307544144 
With 38 sqrt, then 38 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999360966436217 
With 39 sqrt, then 39 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999360966436217 
With 40 sqrt, then 40 times ^2 operations, the number 2.0000000000000000 becomes: 1.9999360966436217 
With 41 sqrt, then 41 times ^2 operations, the number 2.0000000000000000 becomes: 1.9994478907329654 
With 42 sqrt, then 42 times ^2 operations, the number 2.0000000000000000 becomes: 1.9984718365144798 
With 43 sqrt, then 43 times ^2 operations, the number 2.0000000000000000 becomes: 1.9965211562778555 
With 44 sqrt, then 44 times ^2 operations, the number 2.0000000000000000 becomes: 1.9965211562778555 
With 45 sqrt, then 45 times ^2 operations, the number 2.0000000000000000 becomes: 1.9887374575497223 
With 46 sqrt, then 46 times ^2 operations, the number 2.0000000000000000 becomes: 1.9887374575497223 
With 47 sqrt, then 47 times ^2 operations, the number 2.0000000000000000 becomes: 1.9887374575497223 
With 48 sqrt, then 48 times ^2 operations, the number 2.0000000000000000 becomes: 1.9887374575497223 
With 49 sqrt, then 49 times ^2 operations, the number 2.0000000000000000 becomes: 1.8682459487159784 
With 50 sqrt, then 50 times ^2 operations, the number 2.0000000000000000 becomes: 1.6487212645509468 
With 51 sqrt, then 51 times ^2 operations, the number 2.0000000000000000 becomes: 1.6487212645509468 
With 52 sqrt, then 52 times ^2 operations, the number 2.0000000000000000 becomes: 1.0000000000000000 
With 53 sqrt, then 53 times ^2 operations, the number 2.0000000000000000 becomes: 1.0000000000000000 
With 54 sqrt, then 54 times ^2 operations, the number 2.0000000000000000 becomes: 1.0000000000000000 
With 55 sqrt, then 55 times ^2 operations, the number 2.0000000000000000 becomes: 1.0000000000000000 
With 56 sqrt, then 56 times ^2 operations, the number 2.0000000000000000 becomes: 1.0000000000000000 
With 57 sqrt, then 57 times ^2 operations, the number 2.0000000000000000 becomes: 1.0000000000000000 
With 58 sqrt, then 58 times ^2 operations, the number 2.0000000000000000 becomes: 1.0000000000000000 
With 59 sqrt, then 59 times ^2 operations, the number 2.0000000000000000 becomes: 1.0000000000000000 
With 60 sqrt, then 60 times ^2 operations, the number 2.0000000000000000 becomes: 1.0000000000000000 

What is happening is that 1 is returned for n >= 52 as the square root of 2, that is, after 52 times square-root operation, the degree of accuracy required for representing the result goes beyond the degree of accuracy available in a MATLAB double-precision number. Consequently, the later squaring operation on 1.00000000000000 will leave the number unchanged and therefore, 2 is not recovered.

Python

This code (roundoff.py) will yield the following output,

roundoff
With 1 times sqrt and then 1 times **2, the number 2.0 becomes: 2.0000000000000004
With 2 times sqrt and then 2 times **2, the number 2.0 becomes: 1.9999999999999996
With 3 times sqrt and then 3 times **2, the number 2.0 becomes: 1.9999999999999996
With 4 times sqrt and then 4 times **2, the number 2.0 becomes: 1.9999999999999964
With 5 times sqrt and then 5 times **2, the number 2.0 becomes: 1.9999999999999964
With 6 times sqrt and then 6 times **2, the number 2.0 becomes: 1.9999999999999964
With 7 times sqrt and then 7 times **2, the number 2.0 becomes: 1.9999999999999714
With 8 times sqrt and then 8 times **2, the number 2.0 becomes: 2.0000000000000235
With 9 times sqrt and then 9 times **2, the number 2.0 becomes: 2.0000000000000235
With 10 times sqrt and then 10 times **2, the number 2.0 becomes: 2.0000000000000235
With 11 times sqrt and then 11 times **2, the number 2.0 becomes: 2.0000000000000235
With 12 times sqrt and then 12 times **2, the number 2.0 becomes: 1.9999999999991336
With 13 times sqrt and then 13 times **2, the number 2.0 becomes: 1.9999999999973292
With 14 times sqrt and then 14 times **2, the number 2.0 becomes: 1.9999999999973292
With 15 times sqrt and then 15 times **2, the number 2.0 becomes: 1.9999999999973292
With 16 times sqrt and then 16 times **2, the number 2.0 becomes: 2.0000000000117746
With 17 times sqrt and then 17 times **2, the number 2.0 becomes: 2.0000000000408580
With 18 times sqrt and then 18 times **2, the number 2.0 becomes: 2.0000000000408580
With 19 times sqrt and then 19 times **2, the number 2.0 becomes: 2.0000000001573586
With 20 times sqrt and then 20 times **2, the number 2.0 becomes: 2.0000000001573586
With 21 times sqrt and then 21 times **2, the number 2.0 becomes: 2.0000000001573586
With 22 times sqrt and then 22 times **2, the number 2.0 becomes: 2.0000000010885857
With 23 times sqrt and then 23 times **2, the number 2.0 becomes: 2.0000000029511749
With 24 times sqrt and then 24 times **2, the number 2.0 becomes: 2.0000000066771721
With 25 times sqrt and then 25 times **2, the number 2.0 becomes: 2.0000000066771721
With 26 times sqrt and then 26 times **2, the number 2.0 becomes: 1.9999999917775542
With 27 times sqrt and then 27 times **2, the number 2.0 becomes: 1.9999999917775542
With 28 times sqrt and then 28 times **2, the number 2.0 becomes: 1.9999999917775542
With 29 times sqrt and then 29 times **2, the number 2.0 becomes: 1.9999999917775542
With 30 times sqrt and then 30 times **2, the number 2.0 becomes: 1.9999999917775542
With 31 times sqrt and then 31 times **2, the number 2.0 becomes: 1.9999999917775542
With 32 times sqrt and then 32 times **2, the number 2.0 becomes: 1.9999990380770896
With 33 times sqrt and then 33 times **2, the number 2.0 becomes: 1.9999971307544144
With 34 times sqrt and then 34 times **2, the number 2.0 becomes: 1.9999971307544144
With 35 times sqrt and then 35 times **2, the number 2.0 becomes: 1.9999971307544144
With 36 times sqrt and then 36 times **2, the number 2.0 becomes: 1.9999971307544144
With 37 times sqrt and then 37 times **2, the number 2.0 becomes: 1.9999971307544144
With 38 times sqrt and then 38 times **2, the number 2.0 becomes: 1.9999360966436217
With 39 times sqrt and then 39 times **2, the number 2.0 becomes: 1.9999360966436217
With 40 times sqrt and then 40 times **2, the number 2.0 becomes: 1.9999360966436217
With 41 times sqrt and then 41 times **2, the number 2.0 becomes: 1.9994478907329654
With 42 times sqrt and then 42 times **2, the number 2.0 becomes: 1.9984718365144798
With 43 times sqrt and then 43 times **2, the number 2.0 becomes: 1.9965211562778555
With 44 times sqrt and then 44 times **2, the number 2.0 becomes: 1.9965211562778555
With 45 times sqrt and then 45 times **2, the number 2.0 becomes: 1.9887374575497223
With 46 times sqrt and then 46 times **2, the number 2.0 becomes: 1.9887374575497223
With 47 times sqrt and then 47 times **2, the number 2.0 becomes: 1.9887374575497223
With 48 times sqrt and then 48 times **2, the number 2.0 becomes: 1.9887374575497223
With 49 times sqrt and then 49 times **2, the number 2.0 becomes: 1.8682459487159784
With 50 times sqrt and then 50 times **2, the number 2.0 becomes: 1.6487212645509468
With 51 times sqrt and then 51 times **2, the number 2.0 becomes: 1.6487212645509468
With 52 times sqrt and then 52 times **2, the number 2.0 becomes: 1.0000000000000000
With 53 times sqrt and then 53 times **2, the number 2.0 becomes: 1.0000000000000000
With 54 times sqrt and then 54 times **2, the number 2.0 becomes: 1.0000000000000000
With 55 times sqrt and then 55 times **2, the number 2.0 becomes: 1.0000000000000000
With 56 times sqrt and then 56 times **2, the number 2.0 becomes: 1.0000000000000000
With 57 times sqrt and then 57 times **2, the number 2.0 becomes: 1.0000000000000000
With 58 times sqrt and then 58 times **2, the number 2.0 becomes: 1.0000000000000000
With 59 times sqrt and then 59 times **2, the number 2.0 becomes: 1.0000000000000000

What is happening is that 1 is returned for n >= 52 as the square root of 2, that is, after 52 times square-root operation, the degree of accuracy required for representing the result goes beyond the degree of accuracy available in a Python float. Consequently, the later squaring operation on 1.00000000000000 will leave the number unchanged and therefore, 2 is not recovered.

Comments