(phixonline)-->
constant maxIt = 13,
maxItJ = 10
atom a1 = 1.0,
a2 = 0.0,
d1 = 3.2
puts(1," i d\n")
for i=2 to maxIt do
atom a = a1 + (a1 - a2) / d1
for j=1 to maxItJ do
atom x = 0, y = 0
for k=1 to power(2,i) do
y = 1 - 2*y*x
x = a - x*x
end for
a = a - x/y
end for
atom d = (a1-a2)/(a-a1)
printf(1,"%2d %.8f\n",{i,d})
d1 = d
a2 = a1
a1 = a
end for