Greedy algorithm for Egyptian fractions: Difference between revisions

Line 2,202:

The sieve copied from Go

<lang Phix>~~include mpfr.e~~

with javascript_semantics

function egyptian(integer num, denom)

include mpfr.e

mpz n = mpz_init(num),

function egyptian(integer num, denom)

d = mpz_init(denom),

mpz n = mpz_init(num),

t = mpz_init()

d = mpz_init(denom),

sequence result = {}

t = mpz_init()

while mpz_cmp_si(n,0)!=0 do

sequence result = {}

mpz_cdiv_q(t, d, n)

while mpz_cmp_si(n,0)!=0 do

result = append(result,"1/"&mpz_get_str(t))

mpz_cdiv_q(t, d, n)

mpz_neg(d,d)

result = append(result,"1/"&mpz_get_str(t))

mpz_mod(n,d,n)

mpz_neg(d,d)

mpz_neg(d,d)

mpz_mod(n,d,n)

mpz_mul(d,d,t)

mpz_neg(d,d)

end while

mpz_mul(d,d,t)

{n,d} = mpz_free({n,d})

end while

return result

{n,d} = mpz_free({n,d})

end function

return result

end function

procedure efrac(integer num, denom)

string fraction = sprintf("%d/%d",{num,denom}),

procedure efrac(integer num, denom)

prefix = ""

string fraction = sprintf("%d/%d",{num,denom}),

if num>=denom then

prefix = ""

integer whole = floor(num/denom)

if num>=denom then

num -= whole*denom

integer whole = floor(num/denom)

prefix = sprintf("[%d] + ",whole)

num -= whole*denom

end if

prefix = sprintf("[%d] + ",whole)

string e = join(egyptian(num, denom)," + ")

end if

printf(1,"%s -> %s%s\n",{fraction,prefix,e})

string e = join(egyptian(num, denom)," + ")

end procedure

printf(1,"%s -> %s%s\n",{fraction,prefix,e})

end procedure

efrac(43,48)

efrac(5,121)

efrac(43,48)

efrac(2014,59)

efrac(5,121)

efrac(2014,59)

integer maxt = 0,

maxd = 0

integer maxt = 0,

string maxts = "",

maxd = 0

maxds = "",

string maxts = "",

maxda = ""

maxds = "",

maxda = ""

for r=98 to 998 by 900 do -- (iterates just twice!)

sequence sieve = repeat(repeat(false,r+1),r) -- to eliminate duplicates

for r=98 to 998 by 900 do -- (iterates just twice!)

for n=1 to r do

sequence sieve = repeat(repeat(false,r+1),r) -- to eliminate duplicates

for d=n+1 to r+1 do

for n=1 to r do

if sieve[n][d]=false then

for d=n+1 to r+1 do

string term = sprintf("%d/%d",{n,d})

if sieve[n][d]=false then

sequence terms = egyptian(n,d)

string term = sprintf("%d/%d",{n,d})

integer nterms = length(terms)

sequence terms = egyptian(n,d)

if nterms>maxt then

integer nterms = length(terms)

maxt = nterms

if nterms>maxt then

maxts = term

maxt = nterms

elsif nterms=maxt then

~~maxts~~ &= ", " & term

maxts = term

elsif nterms=maxt then

end if

maxts &= ", " & term

integer mlen = length(terms[$])-2

end if

if mlen>maxd then

integer mlen = length(terms[$])-2

maxd = mlen

if mlen>maxd then

maxds = term

maxd = mlen

maxda = terms[$]

maxds = term

elsif mlen=maxd then

maxda = terms[$]

maxds &= ", " & term

elsif mlen=maxd then

end if

maxds &= ", " & term

if n<r/2 then

~~for~~ k=~~2 to 9999~~ do

end if

if n<r/2 then

if d*k > r+1 then exit end if

for k=2 to 9999 do

sieve[n*k][d*k] = true

if d*k > r+1 then exit end if

end for

sieve[n*k][d*k] = true

end if

end for

end if

end if

end for

end if

end for

end for

printf(1,"\nfor proper fractions with 1 to %d digits\n",{length(sprint(r))})

end for

printf(1,"Largest number of terms is %d for %s\n",{maxt,maxts})

printf(1,"\nfor proper fractions with 1 to %d digits\n",{length(sprint(r))})

maxda = maxda[3..$] -- (strip the "1/")

printf(1,"Largest number of terms is %d for %s\n",{maxt,maxts})

maxda[6..-6]="..." -- (show only first/last 5 digits)

maxda = maxda[3..$] -- (strip the "1/")

printf(1,"Largest size for denominator is %d digits (%s) for %s\n",{maxd,maxda,maxds})

maxda[6..-6]="..." -- (show only first/last 5 digits)

end for</lang>

printf(1,"Largest size for denominator is %d digits (%s) for %s\n",{maxd,maxda,maxds})

end for

<pre>