The 147 Puzzle
February 5, 2013
We’ll notate solutions as a list of denominators in ascending order, so the “perfect” solution given above is notated (2,4,7,14,28). The solution is recursive, generalized to k fractions instead of limited to 5:
(define (f k)
(let ((fss (list)))
(let loop1 ((k k) (n 1) (d 1) (fs (list)))
(if (= k 1)
(if (zero? (modulo d n))
(set! fss (cons (cons (/ d n) fs) fss)))
(let ((lo (if (null? fs) (+ (quotient d n) 1)
(max (car fs) (+ (quotient d n) 1))))
(hi (quotient (* k d) n)))
(let loop2 ((a lo))
(when (<= a hi)
(let* ((new-d (lcm a d))
(new-n (- (* n new-d (/ d)) (/ new-d a))))
(loop1 (- k 1) new-n new-d (cons a fs)))
(loop2 (+ a 1)))))))
(reverse (map reverse fss))))
That works, but looks a little bit ugly to me; perhaps someone can suggest a better solution. We keep the numerator and denominator of the accumulating fraction separately, in the variables n and d, because it is easier to compute their quotient when needed than to split the pieces of the fraction; k is the number of fractions remaining, fs is the current list of fractions, and fss is the solutions list. Here’s an example:
> (f 4)
((2 3 7 42) (2 3 8 24) (2 3 9 18) (2 3 10 15) (2 3 12 12)
(2 4 5 20) (2 4 6 12) (2 4 8 8) (2 5 5 10) (2 6 6 6)
(3 3 4 12) (3 3 6 6) (3 4 4 6) (4 4 4 4))
You can run the code at http://programmingpraxis.codepad.org/zxd9grfB, where you will see the 147 solutions to the 5-fraction puzzle. The number of solutions grows very quickly with k: there are 3462 solutions when k is 6, 294314 solutions when k is 7, and 159330691 solutions when k is 8 (A002966). Sums of unit fractions with numerator 1 were used computationally by the ancient Egyptians, so such fractions are sometimes called Egyptian fractions.
Programming Praxis 
[…] today’s Programming Praxis exercise, our goal is to find all combinations of five numbers whose inverses […]
My Haskell solution (see http://bonsaicode.wordpress.com/2013/02/05/programming-praxis-the-147-puzzle/ for a version with comments):
import Data.Ratio puzzle :: Integer -> [[Integer]] puzzle k = f 2 (k%1) 1 where f _ 1 r = if numerator r == 1 then [[denominator r]] else [] f p n r = concat [map (x : ) $ f x (n-1) r' | x <- [p..floor $ n/r], let r' = r - 1%x, r' > 0]A recursive Python version.
from fractions import Fraction import itertools def g(n, low, target, prev): if n == 1: i = 1 / target if i.denominator == 1: print prev + [i.numerator] else: for i in itertools.count(low): fi = Fraction(1, i) if n * fi < target: break t = target - fi if t > 0: g(n-1, i, t, prev + [i]) g(4, 2, 1, [])A cleaner and faster Python version.
from fractions import Fraction def g(n, low, target, prev): invt = 1 / target if n == 1: if invt.denominator == 1: print prev + [invt.numerator] else: for i in range(max(low, int(invt) + 1), int(n * invt) + 1): g(n-1, i, target - Fraction(1, i), prev + [i]) g(4, 2, 1, [])[…] Pages: 1 2 […]
Here’s my Racket take:
http://blog.jverkamp.com/2013/02/06/the-147-puzzle/
Essentially, it ended up being a translation of Paul’s second solution into Racket. My original code was close, but I was missing a bit on the recursive bounds.
The procedure below passes the denominators of each solution to a
procedure which can be used to count or print them, or it can be
carefully crafted to generate the solutions on demand.
I got all confused about the limits and couldn’t work them out until I
took the hint from jpverkamp and copied the limits from Paul. Then my
procedure began to produce sensible results. I’m still not sure why my
lower limit differs by one from his, but it seems to work this way.
Perhaps it’s my obsession with the empty sum.
;;; The procedure calls proc on each different sequence of n ;;; increasing denominators d1, ..., dn dk where dk >= min-d ;;; and (+ (/ d1) ... (/ dn)) = x. (define (for-each-egyptian-sum/minimum-denominator proc x n min-d) (if (and (zero? x) (zero? n)) (proc) (if (positive? x) (do ((dk (max min-d (floor (/ x))) (+ dk 1))) ((< (floor (/ n x)) dk)) (for-each-egyptian-sum/minimum-denominator (lambda rest (apply proc dk rest)) (- x (/ dk)) (- n 1) dk))))) (define (test n) ;; counts the ways 1 can be written as n reciprocals (let ((count 0)) (for-each-egyptian-sum/minimum-denominator (lambda ds (set! count (+ count 1))) 1 n 1) count))I guess recursion is the trickier way to solve this. Why not just iterate?
Also, if you multiply every term by abcde then you get an equation that doesn’t involve division or floats.
for(int a=1; a<10; a++)
{
for(int b=1; b<10; b++)
{
int ab = a*b;
for(int c=1; c<10; c++)
{
int abc = ab*c;
for(int d=1; d<10; d++)
{
int abcd = abc * d;
for(int e=1; e<10; e++)
{
if((b*c*d*e) + (a*c*d*e) + (ab*d*e) + (abc*e) + abcd == abcd*e)
printf("a=%d, b=%d, c=%d, d=%d, e=%d\n", a, b, c, d, e);
}
}
}
}
}
Re why not just iterate: The recursive solutions have the number of terms as a parameter. The same procedure, say f, solves f(1), f(2), f(3), f(4), f(5), and f(6). It solves also more, but after 6 the number of solutions becomes large.
Also, the solution to f(5) containing the largest denominator is (2, 3, 7, 43, 1806), so limiting the range of each denominator to single digits hides solutions. I don’t know if there is a way to know upper bounds to all denominators in advance.
None of the posted solutions use floats. They use exact rational arithmetic.
I think the largest denominator can be calculated (see the function lowest_series below. As you can see the largest denominator becomes very large. In that respect a recursive method works a lot nicer, as the lower and upper bounds can be dermined easily for the remaining terms.
from operator import mul def lowest_series(n): aas = [2, 3] for i in range(3, n + 1): an = reduce(mul, aas) if i < n: an += 1 aas.append(an) return aas for n in range(3, 8): print n, lowest_series(n) #3 [2, 3, 6] #4 [2, 3, 7, 42] #5 [2, 3, 7, 43, 1806] #6 [2, 3, 7, 43, 1807, 3263442] #7 [2, 3, 7, 43, 1807, 3263443, 10650056950806L]