Exercise 1.37 explains that an infinite continued fraction is an expression of the form
The infinite continued fraction where all Ni and Di terms are equal to 1 produces 1/ϕ, where ϕ is the golden ratio described in SICP section 1.2.2.We can approximate the value of an infinite continued fraction by truncating to a given number of terms.
We are asked to define a procedure cont-frac that takes the arguments n, d, and k. n and d are procedures of one argument (the term index i) that compute the Ni and Di terms of the continued fraction. k is the number of terms to expand.We can check our procedure by approximating 1/ϕ using
(cont-frac (lambda (i) 1.0)
(lambda (i) 1.0)
k)
for successive values of
k. We're asked to show how large k must be in order to approximate 1/ϕ to 4 decimal places. As usual, we're also asked to write two versions of the procedure, one that generates a recursive process and one iterative.Since we're given the ending value of
k and it won't do to count down from k to 1 as we compute the accumulated value, we'll define a helper function frac that counts up from 1 to k. I'll show the recursive version first since that comes more naturally.(define (cont-frac n d k)
(define (frac i)
(if (< i k)
(/ (n i) (+ (d i) (frac (+ i 1))))
(/ (n i) (d i))))
(frac 1))
The
cont-frac procedure simply defines the helper function frac then calls it with its starting value of 1. The frac procedure checks to see if i is less than k. If so, it divides the result of applying the n procedure to i by the sum of the result of applying the d procedure to i and a recursive call to frac with the next value of i. If i is equal to k, the recursion ends and one last division operation is performed.We can check the results by using the supplied code with an arbitrary value of
k.> (cont-frac (lambda (i) 1.0)
(lambda (i) 1.0)
5)
0.625
This gives us a value of 1/ϕ that's in the same ball park as the expected value of about 0.61803, but we were asked to find a value of
k where the result is accurate to 4 decimal places. We can just test a few more values to find the answer.> (cont-frac (lambda (i) 1.0)
(lambda (i) 1.0)
8)
0.6176470588235294
> (cont-frac (lambda (i) 1.0)
(lambda (i) 1.0)
9)
0.6181818181818182
> (cont-frac (lambda (i) 1.0)
(lambda (i) 1.0)
10)
0.6179775280898876
We can define an iterative version of
cont-frac by modifying the helper function. Since the iterative procedure will be carrying the intermediate result at each step as a parameter, this time we'll count down from k to 0 instead of counting up from 1 to k as we did before.(define (cont-frac-iter n d k)
(define (frac-iter i result)
(if (= i 0)
result
(frac-iter (- i 1) (/ (n i) (+ (d i) result)))))
(frac-iter (- k 1) (/ (n k) (d k))))
You can run
cont-frac-iter with the same inputs that we used above to verify that you get the same results.Exercise 1.38 asks us to write a program that uses our
cont-frac procedure to approximate e, the base of the natural logarithm. We're to use the continued fraction expansion for e - 2 published by Euler in 1737. In this fraction all of the Ni terms are 1, and the Di terms are1, 2, 1, 1, 4, 1, 1, 6, 1, 1, 8,...
Since
cont-frac was already written in the last exercise, this problem is reduced to writing a good function for computing the Di terms. The series is extremely regular, except for having only one "1" at the beginning. The indices of the non-1 values in the series are 2, 5, 8, 11, 14, 17,... These values are all one less than a multiple of three, or 3i - 1. So we know to begin with that the procedure (d i) should return a 1 when (i + 1) is not divisible by 3.Since the series is so regular, it's easy to find a formula for the remaining values as well. You could do it by plotting the indices and their values in a spreadsheet, but that's not really necessary. First, take a look at the two sequences side-by-side.
i = 2, 5, 8, 11, 14, 17, ...
d(i) = 2, 4, 6, 8, 10, 12, ...
d(i) = 2, 4, 6, 8, 10, 12, ...
What I notice right away is that when I add 3 to the index, the result only increases by 2. This means that I should be able to get the result by dividing the index by 3 and multiplying by 2 (after applying an offset, which is 1 in this case). In mathematical terms:
d(i) = 2(i + 1) / 3
Putting it all together in code, it looks like the following:
(define (d i)
(if (not (= 0 (remainder (+ i 1) 3)))
1
(* 2 (/ (+ i 1) 3))))
Now to compute e, we can use 1 for each Ni term and our new
d procedure for each Di term in a call to cont-frac. Remember that Euler's continued fraction computed e - 2, so we need to take that into account.(define e
(+ 2 (cont-frac (lambda (i) 1.0) d 10)))
> e
2.7182817182817183
Exercise 1.39 asks us to define a procedure
(tan-cf x k) that computes an approximation to the tangent function based on the following continued fraction (published in 1770 by Johann Heinrich Lambert).
The rule for generating the Ni term is fairly simple. If i = 1, the term is equal to x. Otherwise it's -x2. We have to negate all but the first term because our cont-frac procedure adds each term and we need to subtract them in this case.The rule for the Di term is even simpler, since they're just the odd numbers. Di = 2(i) - 1.
(define (square x) (* x x))
(define (tan-cf x k)
(define (n k)
(if (= k 1)
x
(- (square x))))
(define (d k)
(- (* 2 k) 1))
(cont-frac n d k))
We can use Scheme's built-in
tan function to check our work with a few common angles (remember the angle is in radians). We'll stick with k = 10 terms of the continued fraction since we had such good results with that before.> (tan (/ pi 6))
0.5773502691896257
> (tan-cf (/ pi 6) 10)
0.5773502691896257
> (tan (/ pi 4))
0.9999999999999999
> (tan-cf (/ pi 4) 10)
1.0
> (tan (/ pi 3))
1.7320508075688767
> (tan-cf (/ pi 3) 10)
1.732050807568877
Related:
For links to all of the SICP lecture notes and exercises that I've done so far, see The SICP Challenge.
9 comments:
Your code for cont-frac and the procedure to approximate e can be simplified a little.
; Recursive
(define (cont-frac-r n d k)
(define (frac i)
(if (> i k)
0
(/ (n i) (+ (d i) (frac (+ i 1))))))
(frac 1))
; Iterative
(define (cont-frac n d k)
(define (iter i result)
(if (= i 0)
result
(iter (- i 1) (/ (n i) (+ (d i) result)))))
(iter k 0))
; Approximate e
(define (approx-e k)
(define (n i) 1.0)
(define (d i)
(if (not (= 2 (remainder i 3)))
1
(* 2 (+ 1 (quotient i 3)))))
(+ 2 (cont-frac n d k)))
Hi Bill,
Your post has been an inspiration and very helpful in solving SICP. Looking at the answer for 1.38, I find it easier to think in terms of mod, rather than the 2n+1 formula suggested by you.
If there was a 0th element, then every 3rd element is incrementally increasing multiple of 2. Anonymous has also suggested the same. Although adding that we are thinking in terms of mod is what I want to add. This is the code to find value of e - 2.
(cont-frac
(lambda (i) 1.0)
(lambda (i)
(if (= (remainder i 3) 2)
(* 2 (+ (quotient i 3) 1))
1.0))
10)
(Sorry there was a bug in the iterative version before)
Hi, my take at cont-frac is the following (it's Clojure). A little different from all the versions here, but equivalent.
(defn cont-frac [n d k]
(if (= k 0)
(/ (n k) (d k))
(/ (n k) (+ (d k) (cont-frac n d (dec k))))))
(defn cont-frac* [n d k]
(defn iter [k result]
(if (= k 0)
result
(iter (dec k) (/ (n k) (+ (d k) result)))))
(iter (dec k) (/ (n k) (d k))))
Self-correcting myself, my recursive process is not equivalent to yours.
In particular it fails epically with the following, aka when the index order counts:
(defn tan-cf [x k]
(defn n [i] (if (= i 1) x
(- (*' x x))))
(defn d [i] (- (* 2 i) 1))
(cont-frac n d k))
Wanted to share another recursive solution to cont-frac that I found a bit more intuitive:
(define (cont-frac n d k)
(if (= k 0) 0
(/ (n 1) (+ (d 1) (cont-frac
(lambda (i) (n (+ 1 i)))
(lambda (i) (d (+ 1 i)))
(- k 1))))))
The idea is to call yourself with a modified n and d, which represent the "rest of the sequence." I find this more intuitive because it allows you to start at the top of the continued fraction, with N1 and D1, and work your way down recursively. Otherwise you have to start at the bottom.
Bill,
I believe the answer to 1.37 is 11, not 10. As you show, 10 iterations gives 0.6179775280898876 while 11 gives 0.6180555555555556. The question asks for how large you need to make k to get an approximation accurate to 4 decimal places. Since the actual number is ~0.61803402264934 the answer must be 0.6180555555555556.
Thanks.
Post a Comment