(phixonline)-->
with javascript_semantics
constant max_n = 200, branch = 4
atom t0 = time()
include mpfr.e
sequence ivals = repeat(1,2)&repeat(0,max_n-1),
rooted = mpz_inits(max_n+1,ivals),
unrooted = mpz_inits(max_n+1,ivals),
c = mpz_inits(branch,0)
mpz tmp = mpz_init()
procedure tree(integer br, n, l, s, mpz cnt)
mpz cbr = c[br]
for b=br to branch do
s += n
if s>max_n
or (l*2>=s and b>=branch) then
return
end if
if b=br then
mpz_mul(cbr, rooted[n+1], cnt)
else
mpz_add_ui(tmp, rooted[n+1], b-br)
mpz_mul(cbr, cbr, tmp)
mpz_divexact_ui(cbr, cbr, b-br+1)
end if
if l*2<s then
mpz u = unrooted[s+1]
mpz_add(u, u, cbr)
end if
if b<branch then
mpz r = rooted[s+1]
mpz_add(r, r, cbr)
for m=n-1 to 1 by -1 do
tree(b+1, m, l, s, cbr)
end for
end if
end for
end procedure
procedure bicenter(integer s)
if even(s) then
mpz aux = rooted[s/2+1],
u = unrooted[s+1]
mpz_add_ui(tmp, aux, 1)
mpz_mul(tmp, aux, tmp)
mpz_tdiv_q_2exp(tmp, tmp, 1)
mpz_add(u, u, tmp)
end if
end procedure
mpz cnt = mpz_init(1)
integer s = 1
for n=1 to max_n do
tree(1, n, n, s, cnt)
bicenter(n)
if n<=25 or remainder(n,100)=0 then
printf(1,"%d: %s\n",{n, mpz_get_short_str(unrooted[n+1])})
end if
end for
?elapsed(time()-t0)