Methods Based on the Binary Number System
Methods Based on the Binary Number System
The Binary Method
An Egyptian fraction representation can be formed from the binary representation of a number; e.g. 27/22 = 1.0(0111010001)*. (The star indicates that the parenthesized digits repeat in blocks of ten after these initial bits.) The initial nonrepeating part gives fractions 1/2^a, and the repeating part gives fractions 1/(2^a (2^b - 1)), where b is the length of the repetition. We must take some care that the a in the second type of fraction is nonnegative, so we modify the representation above so that there are as many nonrepeating terms as repeating terms: 27/22 = 1.001110100(0101110100)*.
A similar technique works for some other bases than binary. For instance the only digit that causes any trouble in a base 6 representation is 5, but 5/6 = 3/6+2/6 so we can still use this method with base 6. On the other hand this method does not work well with decimal notation as we can not represent 4, 7, 8, or 9 as sums of distinct divisors of 10, so numbers with those digits in their decimal representation would cause a problem for a decimal version of this method.
To implement the binary method, we first define a function to find the binary (or other base) representation of q, returned as two lists of digits. The first member of the first list is the integer part of q, the rest of the first list is the nonrepeating part of the representation, and the second list is the repeating part. It turns out to be easier to find, instead of the digits themselves, certain values mod y from which the digits can be computed. This makes it easier to detect repeating blocks of digits, since they occur exactly when the same value mod y arises twice.
RationalDigits[q_Integer, base_] := {{q},{0}};
RationalDigits[Rational[a_,b_], base_Integer] :=
Module[{nextunit,addunit,units,
reppos,breakpt,finddigit,digitize},
nextunit = (Mod[base Last[#], b]&);
addunit = (Module[{c=nextunit[#]},
If[MemberQ[#,c],
#, Append[#,c]]]&);
units = FixedPoint[addunit, {Mod[a,b]}];
reppos = Position[units, nextunit[units]];
breakpt = reppos[[1]][[1]]-1;
finddigit = (Floor[base # / b]&);
digitize = (finddigit /@ # &);
{Prepend[digitize[Take[units,breakpt]],Floor[a/b]],
digitize[Drop[units,breakpt]]}];
Once we have found the repeating binary representation of a fraction, it is simple to turn the nonzero digits of the representation into terms in an Egyptian fraction representation. Most of the complication in our implementation is due to the point noted earlier, that we should have at least as many nonrepeating digits as repeating digits.
EgyptBinary[q_Integer] := {q};
EgyptBinary[q_Rational] :=
Module[{l = RationalDigits[q,2],
tpow = ({2 #1[[1]], #2}&),
invprod = (#[[2]]/#[[1]]&),
tplist,invlist,
firstlen,firstlist,firstpart,
mul,seclist,secpart,full},
tplist = (FoldList[tpow, {1,#[[1]]}, Drop[#,1]]&);
invlist = (invprod /@ tplist[#])&;
firstlen = Max[Length[l[[1]]],Length[l[[2]]]];
firstlist = Take[Apply[Join,l],firstlen];
firstpart = invlist[firstlist];
mul = 2^Length[l[[2]]]-1;
seclist = RotateRight[l[[2]], Length[l[[1]]]];
secpart = (#/mul& /@ invlist[seclist]);
full = (If[#==0,{ }, #]& /@
Join[firstpart,secpart]);
Flatten[full]];
The correctness of this algorithm is straightforward. The fact that it halts for input q=x/y follows from the fact that the list units computed in RationalDigits is a list of distinct values modulo y, so can have length at most y. The lists of binary digits for x/y have between them at most y elements, and the padding to make the repetition start far enough along at most doubles this. At most y/2 elements on the list can correspond to binary ones, so the eventual Egyption fraction has at most y terms. All denominators are at most 2^(2y).
EgyptBinary[27/22]
1 1 1 1 1 1 1 1 1
{1, -, --, --, ---, ----, ----, -----, -----, ------}
8 16 32 128 2046 8184 16368 32736 130944
The Binary Remainder Method
EBRPower[x_,y_,b_] :=
If[(x*b)~Mod~y < 2*b, b,EBRPower[x,y,2*b]];
EgyptBinRem[Rational[x_,y_]] :=
Module[{p, r, s},
p = EBRPower[x,y,2];
r = Mod[x p, y];
s = Quotient[x p, y];
Join[If[s==0,{},EgyptBinary[s/p]],
If[r==0,{},(#/y&) /@ EgyptBinary[r/p]]]]
The method produces at most Log x + Log y terms; in practice it will typically produce half that many. Each denominator is at most y^2.
EgyptBinRem[18/23]
1 1 1 1
{-, -, --, --}
2 4 46 92
EgyptBinRem[31/311]
1 1 1 1 1 1 1
{--, --, ---, ---, ----, ----, ----}
16 32 311 622 1244 4976 9952
The binary remainder method appears in a paper of Stewart [Ste54], where he uses it to find representations with all denominators even. Similar methods that replace the term p=2^k by some other value have proven useful in many recent results about Egyptian fractions. Breusch [Bre54] and Stewart [Ste54] set p to small multiples of 3^k, to show that every odd-denominator fraction has a representation with all denominators odd. Tenenbaum and Yokota [TY90] use factorial values of p to find representations with (1+o(1))(Log y) / (Log Log y) terms having all denominators bounded by O(y (Log y)^2 / (Log Log y)). Vose [Vos85] uses an even more complicated value of p to show that any x/y has a representation with only O(Sqrt Log[y]) terms. In the generalized remainder method, discussed later, we show how to apply some of these ideas to find short representations.
Egyptian Fractions, Number Theory, David Eppstein, ICS, UC Irvine