Egyptian Fractions
June 4, 2013
One topic that has always fascinated me is the mathematical sophistication of the ancients. The Pythagorean theorem is twenty-five hundred years old. Euclid’s Elements of 300BC remained in use as a geometry textbook for over two millenia, and is still in print today (you can buy the Dover edition for $1). Eratosthenes was sieving prime numbers two centuries before Christ.
Today’s topic is Egyptian fractions, which were mentioned in the Rhind Papyrus of 500BC as an “ancient method” and are thought to date to the time of the pyramids; they survived in use until the 17th or 18th century. Leonardo de Pisa’s famous book Liber Abaci of 1215AD, which introduced to European mathematicians the concept of zero, the Hindu-Arabic system of numeration, and the famous rabbit sequence (Leonardo’s nickname was Fibonacci).
An Egyptian fraction was written as a sum of unit fractions, meaning the numerator is always 1; further, no two denominators can be the same. As easy way to create an Egyptian fraction is to repeatedly take the largest unit fraction that will fit, subtract to find what remains, and repeat until the remainder is a unit fraction. For instance, 7 divided by 15 is less than 1/2 but more than 1/3, so the first unit fraction is 1/3 and the first remainder is 2/15. Then 2/15 is less than 1/7 but more than 1/8, so the second unit fraction is 1/8 and the second remainder is 1/120. That’s in unit form, so we are finished: 7 ÷ 15 = 1/3 + 1/8 + 1/120. There are other algorithms for finding Egyptian fractions, but there is no algorithm that guarantees a maximum number of terms or a minimum largest denominator; for instance, the greedy algorithm leads to 5 ÷ 121 = 1/25 + 1/757 + 1/763309 + 1/873960180913 + 1/1527612795642093418846225, but a simpler rendering of the same number is 1/33 + 1/121 + 1/363.
Your task is to write a program that calculates the ratio of two numbers as an Egyptian fraction. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.
Programming Praxis 
[…] today’s Programming Praxis exercise, our goal is to convert a given fraction to a sum of fractions with […]
My Haskell solution (see http://bonsaicode.wordpress.com/2013/06/04/programming-praxis-egyptian-fractions/ for a version with comments):
import Data.Ratio egypt :: Integer -> Integer -> [Integer] egypt 1 d = [d] egypt n d = e : egypt (numerator r) (denominator r) where (e,r) = (div (d+n-1) n, n%d - 1%e)here’s my general idea.
while numerator > 1:
subtract the largest factor that’s smaller than the numerator; distribute.
if no such factor:
multiply top and bottom by the smallest value that will put the numerator greater than the smallest factor bigger than it.
for example:
7/15 = (2+5)/15 = 1/3 +2/15 = 1/3 +4/30 = 1/3 +1/10 +1/30
25/37 = 50/(2*37) = (37+13)/(2*37) = 1/2 +13/(2*37) = 1/2 +1/37 +11/(2*37) = 1/2 +1/37 +44/(2*4*37) = 1/2 +1/37 +1/8 +7/(2*4*37)
= 1/2 +1/37 +1/8 +1/74 +3/(2*4*37) = 1/2 +1/37 +1/8 +1/74 +1/148 +1/296.
http://pastebin.com/y9KgqXs1
phil’s program loops indefinitely with inputs of 2/4 and 6/8 (perhaps others.) Maybe it needs to reduce fraction to lowest terms before proceeding.
okay thank you david. i tried to fix the program but there seems to be too many road blocks preventing this algorithm form working. currently
31/466 can’t be done.
Clojure solution for a pairing algorithm. Represent x/y as a list of x copies of the fraction 1/y. Reduce the list each step by looking for pairs with the same denominator. Replace each pair with 2/y (when y is even) or 2/(y+1) + 2/(y(y+1)) when y is odd. In each case, the 2 can be factored out, so the numerator will still be 1. The algorithm halts when there are no more equal denominators in the list.
(defn reduce-pairs [l] (cond (< (count l) 2) l (= (first l) (second l)) (let [den (denominator (first l))] (if (even? den) (cons (/ 2 den) (reduce-pairs (rest (rest l)))) (cons (/ 2 (inc den)) (cons (/ 2 (* den (inc den))) (reduce-pairs (rest (rest l))))))) :else (cons (first l) (reduce-pairs (rest l))))) (defn egypt [x y] (loop [l (repeat x (/ 1 y))] (let [l' (sort > (reduce-pairs l))] (if (= l l') l (recur l'))))) user=> (egypt 18 23) (1/2 1/6 1/12 1/35 1/276 1/2415)The Egyptians probably used the pair reduction method because they have a lot of ways in which they would reduce a fraction 2/n where n is odd. but they employed various heuristics to keep the denominators small when n is odd. So this version is also a pair reduction, but uses a lookup table for small odd denominators that looks up the fractions that were used on the Rhind Papyrus. For denominators greater than 101, the identity 1/n + 1/n = 2/(n+1) + 2/(n(n+1) is used.
(def rhind-papyrus [[1/2 1/6], [1/3 1/15], [1/4 1/28] ; 2/3, 2/5, 2/7 [1/6 1/18], [1/6 1/66], [1/8 1/52 1/104] ; 2/9, 2/11, 2/13 [1/10 1/30], [1/12 1/51 1/68], [1/12 1/76 1/114] ; 2/15, 2/17, 2/19 [1/14 1/42], [1/12 1/276], [1/15 1/75] ; 2/21, 2/23, 2/25 [1/18 1/54], [1/24 1/58 1/174 1/232] ; 2/27 2/29 [1/20 1/124 1/155], [1/22 1/66], [1/30 1/42] ; 2/31, 2/33, 2/35 [1/24 1/111 1/296], [1/26 1/78], [1/24 1/246 1/328] ; 2/37, 2/39, 2/41 [1/42 1/86 1/129 1/301], [1/30 1/90] ; 2/43, 2/45 [1/30 1/141 1/470], [1/28 1/196], [1/34 1/102] ; 2/47, 2/49, 2/51 [1/30 1/318 1/795], [1/30 1/330], [1/38 1/114] ; 2/53, 2/55, 2/57 [1/36 1/236 1/531], [1/40 1/244 1/488 1/610] ; 2/59, 2/61 [1/42 1/126], [1/39 1/195], [1/40 1/335 1/536] ; 2/63, 2/65, 2/67 [1/46 1/138], [1/40 1/568 1/710], [1/60 1/219 1/292 1/365] ; 2/69, 2/71 2/73 [1/50 1/150], [1/44 1/308], [1/60 1/237 1/316 1/790] ; 2/75, 2/77, 2/79 [1/54 1/162], [1/60 1/332 1/415 1/498] ; 2/81, 2/83 [1/39 1/195], [1/58 1/174], [1/60 1/356 1/534 1/890] ; 2/85, 2/87, 2/89 [1/70 1/130], [1/62 1/186], [1/60 1/380 1/570] ; 2/91 2/93 2/95 [1/56 1/679 1/776], [1/66 1/198], [1/101 1/202 1/303 1/606]] ; 2/97 2/99 2/101 ) (defn ahmes [n] (rhind-papyrus (dec (/ (dec n) 2)))) (defn reduce-pairs [l] (cond (< (count l) 2) l (= (first l) (second l)) (let [den (denominator (first l))] (cond (even? den) (cons (/ 2 den) (reduce-pairs (rest (rest l)))) (<= den 101) (reduce conj (ahmes den) (reduce-pairs (rest (rest l)))) :else (cons (/ 2 (inc den)) (cons (/ 2 (* den (inc den))) (reduce-pairs (rest (rest l))))))) :else (cons (first l) (reduce-pairs (rest l))))) (defn egypt [x y] (loop [l (repeat x (/ 1 y))] (let [l' (sort > (reduce-pairs l))] (if (= l l') l (recur l'))))) user=> (egypt 18 23) (1/2 1/6 1/12 1/46 1/138 1/276)There is a binary method which takes advantage of the fact that a non-repeating fraction in binary notation is already a sum of unit fractions (e.g. 5/16 = 0.0101 = 1/16 + 1/4, so the method splits a fraction x/y into a sum of inverse powers of two and a remainder portion that is also inverse powers of two multiplied by 1/y. See http://www.ics.uci.edu/~eppstein/numth/egypt/binary.html. Cook thing is that the number of terms is bounded by log(x) + log(y) and denominator is no larger than y^2.
(def powers_of_two (iterate (partial * 2) 1)) (defn bin_frac [x y] ; create Egyptian fraction from x/y; y MUST BE power of 2 (loop [x x, y y, result ()] (cond (= x 0) result (odd? x) (recur (quot x 2) (quot y 2) (cons (/ 1 y) result)) :else (recur (quot x 2) (quot y 2) result)))) (defn bin_divisor [x y] (first (drop-while (fn [p] (>= (mod (* x p) y) (* p 2))) powers_of_two))) (defn egypt [x y] (if (> x y) (cons (quot x y) (egypt (rem x y) y)) ;else (let [p (bin_divisor x y), r (rem (* x p) y), s (quot (* x p) y)] (concat (bin_frac s p) (map #(/ %1 y) (bin_frac r p)))))) user=> (egypt 5 121) (1/32 1/121 1/968 1/1936 1/3872) user=> (egypt 355 113) ; pi (3 1/8 1/113 1/226 1/452 1/904) user=> (egypt 18 23) (1/2 1/4 1/46 1/92)[…] This post presents C# solutions to a coin change problem as described in https://programmingpraxis.com/2013/06/04/egyptian-fractions. […]