SICP-2nd-Edition-Exercise-1.7

From Boozled

Jump to: navigation, search

Contents

Exercise 1.7 from SICP

The good-enough? test used in computing square roots will not be very effective for finding the square roots of very small numbers. Also, in real computers, arithmetic operations are almost always performed with limited precision. This makes our test inadequate for very large numbers. Explain these statements, with examples showing how the test fails for small and large numbers. An alternative strategy for implementing good-enough? is to watch how guess changes from one iteration to the next and to stop when the change is a very small fraction of the guess. Design a square-root procedure that uses this kind of end test. Does this work better for small and large numbers?

The question is using the Babylonian Method to compute square roots.

x_0 \approx \sqrt{S}
x_{n+1} = \frac{1}{2} \left(x_n + \frac{S}{x_n}\right)
\sqrt S = \lim_{n \to \infty} x_n

Some notation to be aware of:

ε = acceptable error

Attempt

Small Numbers

I am not actually too sure what the authors wanted here. I don't expect they wanted us to try and calculate a root where ε is larger than our root. I have assumed they want us to compute ever lower numbers until we hit a problem then explain it.

If we make the number we are seeking the root for very small we must also make the ε small enough to be able to deal with it. For instance if we use 

  1. (define v 0.0000001)

there we must also change the definition of our good-enough? procedure accordingly. For instance the following procedure will compute a bad result if ε is not sufficiently smaller than the number we are seeking root to. 

  1. (define (good-enough? guess x)
  2.   (define tmp "")
  3.   (set! tmp (abs (- (square guess) x))) ; This extra step is here for debugging!
  4.   (< tmp epsilon))

So we need a procedure to make sure that ε is sufficiently small enough that we can compute the result to see if there is a problem as we decrease the root sought. 

  1. (define (good-enough? guess x)
  2.   (define tmp "")
  3.   (set! tmp (abs (- (square guess) x)))
  4.   (< tmp epsilon))
  5.  
  6. (define sought 0.00000000001)
  7. (define epsilon (/ sought 10))
  8. (display (mysqrt sought 1.0))

If I use a standard calculator to give me the answer to sqrt(0.00000000001 I get 3.16227766e-06. The function above gives me 3.170623279538559e-06. Our relative error can be calculated as follows.

\eta = \frac{|v_{\text{approx}}-v|}{|v|},
0.00263807 \approx \frac{|3.170623279538559\times10^-6 - 3.16227766\times10^-6|}{3.16227766\times10^-6},

We can improve on our relative error by changing the value of (epsilon) in our program.

Going back to the point that I cannot see what the authors wanted here because I seem to be able to make the function quite accurate for small values by changing our required ε. I assume they wanted us to point out that the original ε in the book was not small enough to compute accurate roots for values approaching or less than that number. I could be wrong here.

If we decrease the number close to ε we will most certainly get inaccurate results. For instance if ε = 0.0001 is fixed:

X Relative Error
0.001 0.03549633724272905
0.0008 0.05496301527701619
0.0006 0.0037842376707752923
0.0004 0.012015644103529066
0.0002 0.055003946635965306
0.00015 0.090908194076472
0.00013 0.11414024646177157
0.00012 0.12888245342292026
0.00011 0.12888245342292026
0.0001 0.1675473894984787


We can see that when our X = ε we have a relative error of 0.17. This is not particularly great.

Large Numbers

For large numbers we can see we run into problems trying to represent floating point numbers. For instance, given:

ε = 0.01
X = 70368769910000

we got the following output

guess=8394239.560662625 x=70368769910000 epsilon= 0.01 tmp= 94487891793.46875
guess=8388611.42179962 x=70368769910000 epsilon= 0.01 tmp= 31675947.046875
guess=8388609.53376682 x=70368769910000 epsilon= 0.01 tmp= 3.578125
guess=8388609.533766605 x=70368769910000 epsilon= 0.01 tmp= 0.015625
guess=8388609.533766605 x=70368769910000 epsilon= 0.01 tmp= 0.015625

Our value for guess is now repeating so we have entered an infinite loop. For more information describing what is going on there have a look at

Wikipedia Floating Point.
Floating-Point Representation at David Russinoffs site


Alternative Method

The question asks:

An alternative strategy for implementing good-enough? is to watch how guess changes from one iteration to the next and to stop when the change is a very small fraction of the guess.


So we need to track how our guess changes and determine the ratio between iterations. This is simple enough if we resort to using global variables but I would prefer not to use this method, somehow using global variables does not seem like a Scheme/Lisp idiom I should be employing for this. The following code is my attempt at this. 

  1. (define (sqrt-iter2 guess x)
  2.   (define next-guess (improve guess x))
  3.   (if (good-enough2? guess next-guess)
  4.       guess      
  5.       (sqrt-iter2 next-guess x)
  6.       ))
  7.  
  8. (define (good-enough2? last-guess this-guess) 
  9.   (define tmp "")
  10.   (set! tmp (abs(/ this-guess last-guess)))
  11.   (display "this-guess=") (display this-guess) (display " last-guess=") (display last-guess) (display " tmp=") (display tmp)  (newline)
  12.   (> epsilon (abs (- 1 tmp))))

The method above highlights something that may be of interest. If ε = 0.01 then for X = 70368769910000 we get a relative error of err \approx 0.0006711514. This is accurate but when we are using numbers as large as this if we were using this in a chain of calculations it can have a large impact. For instance, drscheme gives me (display (sqrt 70368769910000 )) = 8388609.533766607. Ysing this we can see what this equates in real terms (8388611.42179962 - 8388609.533766607)^2 \approx 3.56. If we had hundreds of calculations happening then this error could have quite an impact.

This new method works with numbers that are both much larger and much smaller than the previous attempts and seems to work well. At least this code completes for all the numbers I have tried with it.

Scheme

The following code may need some tidying up but it contains all the procedure I've used to complete this Exercise in their raw form. Some of the earlier procedures have been commented out etc. If you are cutting and pasting this code then you need to ask if you are actually doing the Exercise or not.

Thought: Perhaps I should only provide code snippets to discourage cut&paste. 

  1. (define (square x) ( * x x))
  2. (define (average x y)
  3.   (/ (+ x y) 2))
  4.  
  5. (define (improve guess x)
  6.   (average guess (/ x guess)))
  7.  
  8. ;(define (sqrt-iter guess x)
  9. ;  (if (good-enough? guess x)
  10. ;      guess      
  11. ;      (sqrt-iter (improve guess x) x)
  12. ;      ))
  13.  
  14. (define (sqrt-iter2 guess x)
  15.   (define next-guess (improve guess x))
  16.   (if (good-enough2? guess next-guess)
  17.       guess      
  18.       (sqrt-iter2 next-guess x)
  19.       ))
  20.  
  21. (define (good-enough2? last-guess this-guess) 
  22.   (define tmp "")
  23.   (set! tmp (abs(/ this-guess last-guess)))
  24.   ; (display "this-guess=") (display this-guess) (display " last-guess=") (display last-guess) (display " tmp=") (display tmp)  (newline)
  25.   (> epsilon (abs (- 1 tmp))))
  26.  
  27. ;(define (good-enough? guess x)
  28. ;  (< (abs (- (square guess) x)) epsilon))
  29. ; The following function has side effects ie the display statements!
  30. ;(define (good-enough? guess x)
  31. ;  (define tmp "")
  32. ;  (set! tmp (abs (- (square guess) x)))
  33. ;  (display "guess=") (display guess) (display " x=") (display x) (display " epsilon= ") (display epsilon)  (display " tmp= ") (display tmp) (newline)
  34. ;  (< tmp epsilon))
  35.  
  36. (define (mysqrt x start)
  37.     (sqrt-iter2 start x))
  38.  
  39.  
  40. (define (abs-error approx actual) 
  41.   (abs (- approx actual)))
  42.  
  43. (define (rel-error approx actual) 
  44.   (- (/ approx actual) 1))
  45.  
  46. ;(define sought 70368744180000)
  47. (define sought 100)
  48. (define epsilon 0.00001);(/ sought 10))
  49.  
  50. (define sq (mysqrt sought 1.0))
  51. (display sq)
  52. (newline)
  53. (define re (rel-error (mysqrt sought 1.0) (sqrt sought)) )
  54. (newline)
  55. (display "Rel Error: ") 
  56. (display re) 
  57. (newline)

Parent Course

6.001_Structure_and_Interpretation_of_Computer_Programs

Further Reading

Computing Roots at Wikipedia
Wikipedia Floating Point.
Floating-Point Representation at David Russinoffs site
Retrieved from "http://www.boozled.org/wiki/index.php/SICP-2nd-Edition-Exercise-1.7"
Personal tools