Zeckendorf arithmetic: Difference between revisions

Line 1,414:
XOXOOOXOX *z XOXOOX = XOXOXOOOOOXOXOX # multiplication
XOXOXOOOOOXOXOX /z XOO = XOOXOXOOOXOOX remainder XO # division
</pre>
 
=={{header|Phix}}==
Uses a binary representation of Zeckendorf numbers, eg decimal 11 is stored as 0b10100, ie meaning 8+3, but actually 20 in decimal.<br>
As such, they can be directly compared using the standard comparison operators, and printed quite trivially just by using the %b format.<br>
They are however (and not all that surprisingly) pulled apart into individual bits for addition/subtraction, etc.<br>
Does not handle negative numbers or anything >139583862445 (-ve probably doable but messy, >1.4e12 requires a total rewrite, probably using string representation).
<lang Phix>sequence fib = {1,1}
 
function zeckendorf(atom n)
-- Same as [[Zeckendorf_number_representation#Phix]]
atom r = 0
while fib[$]<n do
fib &= fib[$] + fib[$-1]
end while
integer k = length(fib)
while k>2 and n<fib[k] do
k -= 1
end while
for i=k to 2 by -1 do
integer c = n>=fib[i]
r += r+c
n -= c*fib[i]
end for
return r
end function
 
function decimal(object z)
-- Convert Zeckendorf number(s) to decimal
atom dec = 0, bit = 2
if sequence(z) then
for i=1 to length(z) do
z[i] = decimal(z[i])
end for
return z
end if
while z do
if and_bits(z,1) then
dec += fib[bit]
end if
bit += 1
if bit>length(fib) then
fib &= fib[$] + fib[$-1]
end if
z = floor(z/2)
end while
return dec
end function
 
function to_bits(integer x)
-- Simplified copy of int_to_bits(), but in reverse order,
-- and +ve only but (also only) as many bits as needed, and
-- ensures there are *two* trailing 0 (most significant)
sequence bits = {}
if x<0 then ?9/0 end if -- sanity/avoid infinite loop
while 1 do
bits &= remainder(x,2)
if x=0 then exit end if
x = floor(x/2)
end while
bits &= 0 -- (since eg 101+101 -> 10000)
return bits
end function
 
function to_bits2(integer a,b)
-- Apply to_bits() to a and b, and pad to the same length
sequence sa = to_bits(a), sb = to_bits(b)
integer diff = length(sa)-length(sb)
if diff!=0 then
if diff<0 then sa &= repeat(0,-diff)
else sb &= repeat(0,+diff)
end if
end if
return {sa,sb}
end function
 
function to_int(sequence bits)
-- Copy of bits_to_int(), but in reverse order (lsb last)
atom val = 0, p = 1
for i=length(bits) to 1 by -1 do
if bits[i] then
val += p
end if
p += p
end for
return val
end function
function zstr(object z)
if sequence(z) then
for i=1 to length(z) do
z[i] = zstr(z[i])
end for
return z
end if
return sprintf("%b",z)
end function
 
function rep(sequence res, integer ds, sequence was, wth)
-- helper for cleanup, validates replacements
integer de = ds+length(was)-1
if res[ds..de]!=was then ?9/0 end if
if length(was)!=length(wth) then ?9/0 end if
res[ds..de] = wth
return res
end function
 
function zcleanup(sequence res)
-- (shared by zadd and zsub)
integer l = length(res)
-- first stage, left to right, {020x -> 100x', 030x -> 110x', 021x->110x, 012x->101x}
for i=1 to l-3 do
switch res[i..i+2]
case {0,2,0}: res[i..i+2] = {1,0,0} res[i+3] += 1
case {0,3,0}: res[i..i+2] = {1,1,0} res[i+3] += 1
case {0,2,1}: res[i..i+2] = {1,1,0}
case {0,1,2}: res[i..i+2] = {1,0,1}
end switch
end for
-- first stage cleanup
if l>1 then
if res[l-1]=3 then res = rep(res,l-2,{0,3,0},{1,1,1}) -- 030 -> 111
elsif res[l-1]=2 then
if res[l-2]=0 then res = rep(res,l-2,{0,2,0},{1,0,1}) -- 020 -> 101
else res = rep(res,l-3,{0,1,2,0},{1,0,1,0}) -- 0120 -> 1010
end if
end if
end if
if res[l]=3 then res = rep(res,l-1,{0,3},{1,1}) -- 03 -> 11
elsif res[l]=2 then
if res[l-1]=0 then res = rep(res,l-1,{0,2},{1,0}) -- 02 -> 10
else res = rep(res,l-2,{0,1,2},{1,0,1}) -- 012 -> 101
end if
end if
-- second stage, pass 1, right to left, 011 -> 100
for i=length(res)-2 to 1 by -1 do
if res[i..i+2]={0,1,1} then res[i..i+2] = {1,0,0} end if
end for
-- second stage, pass 2, left to right, 011 -> 100
for i=1 to length(res)-2 do
if res[i..i+2]={0,1,1} then res[i..i+2] = {1,0,0} end if
end for
return to_int(res)
end function
 
function zadd(integer a, b)
sequence {sa,sb} = to_bits2(a,b)
return zcleanup(reverse(sq_add(sa,sb)))
end function
 
function zinc(integer a)
return zadd(a,0b1)
end function
 
function zsub(integer a, b)
sequence {sa,sb} = to_bits2(a,b)
sequence res = reverse(sq_sub(sa,sb))
-- (/not/ combined with the first pass of the add routine!)
for i=1 to length(res)-2 do
switch res[i..i+2] do
case {1, 0, 0}: res[i..i+2] = {0,1,1}
case {1,-1, 0}: res[i..i+2] = {0,0,1}
case {1,-1, 1}: res[i..i+2] = {0,0,2}
case {1, 0,-1}: res[i..i+2] = {0,1,0}
case {2, 0, 0}: res[i..i+2] = {1,1,1}
case {2,-1, 0}: res[i..i+2] = {1,0,1}
case {2,-1, 1}: res[i..i+2] = {1,0,2}
case {2, 0,-1}: res[i..i+2] = {1,1,0}
end switch
end for
-- copied from PicoLisp: {1,-1} -> {0,1} and {2,-1} -> {1,1}
for i=1 to length(res)-1 do
switch res[i..i+1] do
case {1,-1}: res[i..i+1] = {0,1}
case {2,-1}: res[i..i+1] = {1,1}
end switch
end for
if find(-1,res) then ?9/0 end if -- sanity check
return zcleanup(res)
end function
 
function zdec(integer a)
return zsub(a,0b1)
end function
 
function zmul(integer a, b)
integer res = 0
sequence mult = {a,zadd(a,a)} -- (as per task desc)
integer bits = 2
while bits<b do
mult = append(mult,zadd(mult[$],mult[$-1]))
bits *= 2
end while
integer bit = 1
while b do
if and_bits(b,1) then
res = zadd(res,mult[bit])
end if
b = floor(b/2)
bit += 1
end while
return res
end function
 
function zdiv(integer a, b)
integer res = 0
sequence mult = {b,zadd(b,b)}
integer bits = 2
while mult[$]<a do
mult = append(mult,zadd(mult[$],mult[$-1]))
bits *= 2
end while
for i=length(mult) to 1 by -1 do
integer mi = mult[i]
if mi<=a then
res = zadd(res,bits)
a = zsub(a,mi)
if a=0 then exit end if
end if
bits = floor(bits/2)
end for
return {res,a} -- (a is the remainder)
end function
 
for i=0 to 20 do
integer zi = zeckendorf(i)
atom d = decimal(zi)
printf(1,"%2d: %7b (%d)\n",{i,zi,d})
end for
 
procedure test(atom a, string op, atom b, object res, string expected)
string zres = iff(atom(res)?zstr(res):join(zstr(res)," rem ")),
dres = sprintf(iff(atom(res)?"%d":"%d rem %d"),decimal(res)),
aka = sprintf("aka %d %s %d = %s",{decimal(a),op,decimal(b),dres}),
ok = iff(zres=expected?"":" *** ERROR ***!!")
printf(1,"%s %s %s = %s, %s %s\n",{zstr(a),op,zstr(b),zres,aka,ok})
end procedure
 
test(0b0,"+",0b0,zadd(0b0,0b0),"0")
test(0b101,"+",0b101,zadd(0b101,0b101),"10000")
test(0b10100,"-",0b1000,zsub(0b10100,0b1000),"1001")
test(0b100100,"-",0b1000,zsub(0b100100,0b1000),"10100")
test(0b1001,"*",0b101,zmul(0b1001,0b101),"1000100")
test(0b1000101,"/",0b101,zdiv(0b1000101,0b101),"1001 rem 1")
 
test(0b10,"+",0b10,zadd(0b10,0b10),"101")
test(0b101,"+",0b10,zadd(0b101,0b10),"1001")
test(0b1001,"+",0b1001,zadd(0b1001,0b1001),"10101")
test(0b10101,"+",0b1000,zadd(0b10101,0b1000),"100101")
test(0b100101,"+",0b10101,zadd(0b100101,0b10101),"1010000")
test(0b1000,"-",0b101,zsub(0b1000,0b101),"1")
test(0b10101010,"-",0b1010101,zsub(0b10101010,0b1010101),"1000000")
test(0b1001,"*",0b101,zmul(0b1001,0b101),"1000100")
test(0b101010,"+",0b101,zadd(0b101010,0b101),"1000100")
 
test(0b10100,"+",0b1010,zadd(0b10100,0b1010),"101000")
test(0b101000,"-",0b1010,zsub(0b101000,0b1010),"10100")
 
test(0b100010,"*",0b100101,zmul(0b100010,0b100101),"100001000001")
test(0b100001000001,"/",0b100,zdiv(0b100001000001,0b100),"101010001 rem 0")
test(0b101000101,"*",0b101001,zmul(0b101000101,0b101001),"101010000010101")
test(0b101010000010101,"/",0b100,zdiv(0b101010000010101,0b100),"1001010001001 rem 10")
 
test(0b10100010010100,"+",0b1001000001,zadd(0b10100010010100,0b1001000001),"100000000010101")
test(0b10100010010100,"-",0b1001000001,zsub(0b10100010010100,0b1001000001),"10010001000010")
test(0b10000,"*",0b1001000001,zmul(0b10000,0b1001000001),"10100010010100")
test(0b1010001010000001001,"/",0b100000000100000,zdiv(0b1010001010000001001,0b100000000100000),"10001 rem 10100001010101")
 
test(0b10100,"+",0b1010,zadd(0b10100,0b1010),"101000")
test(0b10100,"-",0b1010,zsub(0b10100,0b1010),"101")
test(0b10100,"*",0b1010,zmul(0b10100,0b1010),"101000001")
test(0b10100,"/",0b1010,zdiv(0b10100,0b1010),"1 rem 101")
integer m = zmul(0b10100,0b1010)
test(m,"/",0b1010,zdiv(m,0b1010),"10100 rem 0")</lang>
{{out}}
<pre>
0: 0 (0)
1: 1 (1)
2: 10 (2)
3: 100 (3)
4: 101 (4)
5: 1000 (5)
6: 1001 (6)
7: 1010 (7)
8: 10000 (8)
9: 10001 (9)
10: 10010 (10)
11: 10100 (11)
12: 10101 (12)
13: 100000 (13)
14: 100001 (14)
15: 100010 (15)
16: 100100 (16)
17: 100101 (17)
18: 101000 (18)
19: 101001 (19)
20: 101010 (20)
0 + 0 = 0, aka 0 + 0 = 0
101 + 101 = 10000, aka 4 + 4 = 8
10100 - 1000 = 1001, aka 11 - 5 = 6
100100 - 1000 = 10100, aka 16 - 5 = 11
1001 * 101 = 1000100, aka 6 * 4 = 24
1000101 / 101 = 1001 rem 1, aka 25 / 4 = 6 rem 1
10 + 10 = 101, aka 2 + 2 = 4
101 + 10 = 1001, aka 4 + 2 = 6
1001 + 1001 = 10101, aka 6 + 6 = 12
10101 + 1000 = 100101, aka 12 + 5 = 17
100101 + 10101 = 1010000, aka 17 + 12 = 29
1000 - 101 = 1, aka 5 - 4 = 1
10101010 - 1010101 = 1000000, aka 54 - 33 = 21
1001 * 101 = 1000100, aka 6 * 4 = 24
101010 + 101 = 1000100, aka 20 + 4 = 24
10100 + 1010 = 101000, aka 11 + 7 = 18
101000 - 1010 = 10100, aka 18 - 7 = 11
100010 * 100101 = 100001000001, aka 15 * 17 = 255
100001000001 / 100 = 101010001 rem 0, aka 255 / 3 = 85 rem 0
101000101 * 101001 = 101010000010101, aka 80 * 19 = 1520
101010000010101 / 100 = 1001010001001 rem 10, aka 1520 / 3 = 506 rem 2
10100010010100 + 1001000001 = 100000000010101, aka 888 + 111 = 999
10100010010100 - 1001000001 = 10010001000010, aka 888 - 111 = 777
10000 * 1001000001 = 10100010010100, aka 8 * 111 = 888
1010001010000001001 / 100000000100000 = 10001 rem 10100001010101, aka 9876 / 1000 = 9 rem 876
10100 + 1010 = 101000, aka 11 + 7 = 18
10100 - 1010 = 101, aka 11 - 7 = 4
10100 * 1010 = 101000001, aka 11 * 7 = 77
10100 / 1010 = 1 rem 101, aka 11 / 7 = 1 rem 4
101000001 / 1010 = 10100 rem 0, aka 77 / 7 = 11 rem 0
</pre>
 
7,806

edits