sicp-ex-1.31

<< Previous exercise (1.30) | sicp-solutions | Next exercise (1.32) >>

a.

   (define (product term a next b) 
     (if (> a b) 1 
        (* (term a) (product term (next a) next b)))) 

factorial function, in terms of product function above, can be written as below.

  
 (define (identity x) x) 
  
 (define (next x) (+ x 1)) 
  
 (define (factorial n) 
 (product identity 1 next n)) 

approximations to pi using john wallis' formula, can be found by defining a new term to be used by product function as below

 (define (pi-term n) 
   (if (even? n) 
       (/ (+ n 2) (+ n 1)) 
       (/ (+ n 1) (+ n 2)))) 

And it can be used as follows.

 (* (product pi-term 1 next 6) 4)   ;;= 3.3436734693877552 
 (* (product pi-term 1 next 100) 4) ;;= 3.1570301764551676 

b.

  (define (product term a next b) 
   (define (iter a res) 
     (if (> a b) res 
          (iter (next a) (* (term a) res)))) 
    (iter a 1)) 

Another approach would be to recognize the series as 4*4*6*6*.../3*3*5*5*... The leading factor of 2 in the numerator must be dealt with.

 (define (pi n) 
   (define (term x) (* x x)) 
   (define (next x) (+ x 2)) 
   ; since we are increasing the numbers by two on every iteration 
   (define limit (* n 2)) 
   ; upper term: - 2 always goes first, start building product from 4 
   ;             - as 2 numbers are skipped, the limit must respect that, too 
   ;             - since we are squaring one time too often at the end, 
   ;               we have to divide that back out of the result 
   ; lower term: start with 3, which is 1 more than the upper term 
   ;             -> so increase limit by 1 
   ;  
   (* 4 (/ (/ (* 2 (product term 4 next (+ limit 2))) 
              (+ limit 2)) 
           (product term 3 next (+ limit 1))))) 

Yielding:

 (pi 6) ;;= 3.255721745 
 (pi 100) ;;= 3.149378473 

n defines steps in terms of pairs of factors. So where the above solution does 6 steps when n=6, this solution does 12 steps when n=6:

 (pi 3) ;;= 3,343673469 
 (pi 50) ;;= 3,157030176 

A third way of defining the pi-term function, using integer arithmetic instead of the even? predicate:

 (define (pi n) 
   (define (pi-term n) 
     (/ (* 2 (ceiling (/ (+ 1 n) 2))) 
        (+ 1 (* 2 (ceiling (/ n 2)))))) 
   (product pi-term 1 inc n)) 
  

A fourth way of doing it is to recognize pi = 4 * (2 * (1/3) *4 * (1/5) ...) * ((1/3) * 4 * (1/5) ...)

Arithmetic progression formula was used to generate the nth number of the sequence

namely a_k = a_1 + (k-1)d

 (define (pi-prod n) 
   (define (inv k) 
     (if (even? k) 
         k 
         (/ 1 k))) 
  (* 4.0 (* (product 2 (+ 2 (- n 1)) inc inv ) 
               (product 3 (+ 3 (- n 1)) inc inv)))) 

Another way is to recognise that the formula for pi is composed of two separate products:

pi = 4 * ((2/3 * 4/5 * 6/7 ... * n/(n+1)) * (4/3 * 6/5 * ... * n/(n-1))

This gives:

 (define (approx-pi n) 
   (* 4 
    (product 
    (lambda (n) (/ n (+ 1 n))) 
    2 
    (lambda (n) (+ n 2)) 
    n) 
    (product 
    (lambda (n) (/ n (- n 1))) 
    4 
    (lambda (n) (+ n 2)) 
    n))) 

 (define (product-iter f a next b) 
   (define (iter a result) 
     (if (> a b) 
         result 
         (iter (next a) (* (f a) result)))) 
   (iter a 1)) 
  
 (define (pi-term n) 
   (* (/ (- n 1.0) n) (/ (+ n 1.0) n))) 
  
 (define (next-pi n) 
   (+ n 2.0)) 
  
 (define (pi accuracy) 
       (* 4 (product-iter pi-term 3.0 next-pi accuracy))) 

My method is the same as the above. A clarification though, this method works by splitting the formula into:

(2⋅4/3⋅3) ⋅ (4⋅6/5⋅5) ⋅ (6⋅8/7⋅7)

We can then see the general expression for each pair of terms:

((n-1)⋅(n+1))/n^2 --> (n^2 - 1)/n^2 --> n^2/n^2 - 1/n^2 --> 1 - 1/n^2

This translates into code as:

 (define (pi-term a) 
     (- 1 (/ 1 (square a)))) 
  
 (define (pi-next a) 
     (+ a 2)) 
  
 (define (pi accuracy) 
     (* 4.0 
        (product-iter pi-term 3 pi-next (+ (* accuracy 2) 1)))) 

The method I used for the Wallis Product recognizes the pairs as following:

(2/3) -> 2.5

(4/3) -> 3.5

(4/5) -> 4.5

(6/5) -> 5.5

...

We first of all define our sum (named Pi as in the mathematical notation):

 (define (pi term next n upper_bound) 
         (if (> n upper_bound) 
             1 
             (* (term n) 
                (pi term next (next n) upper_bound) 
             ) 
         ) 
 ) 

Then we define the Wallis Product:

 (define (wallis a b pi) 
         (define (wallis-term x) 
                 (if (= (remainder (inexact->exact (floor x)) 2) 0) 
                     (/ (floor x) (ceiling x)) 
                     (/ (ceiling x) (floor x)) 
                 ) 
         ) 
          
         (define (wallis-next x) 
                 (+ x 1) 
         ) 
          
         (* 4 (pi wallis-term wallis-next a b)) 
 ) 

And then we can test it:

 (display (wallis 2.5 1000 pi))