Hello Maxima users,
I was in need of continued fraction expansion of square roots, and there is something weird:
cf(sqrt(n)) returns (apparently) what is expected except for some numbers:
?? sublist(1 .. 1000, lambda([n], integerp(cf(sqrt(n)))));
?? [8,27,32,72,108,128,200,243,288,392,432,512,648,675,800,968,972]
Where a .. b is defined by
?? ".."(a, b) := makelist(i, i, a, b, if a < b then 1 else -1)$
?? infix("..")$
That is, for these numbers, the result isn't even a list.
And the result changes with cflength (for cflength > 1, I check numbers 1 .. 10000 to get more answers)
cflength=1:??? [8,27,32,72,108,128,200,243,288,392,432,512,648,675,800,968,972]
cflength=2:??? [200,800,1800,3200,3267,5000,7200,9800]
cflength=3:??? [288,1152,2592,4608,7200]
cflength=4:??? [6728]
cflength=5:??? [9800]
Now, it's not too important in my case, because the algorithm is quite simple for sqrt(n):
cf_sqrt(d) := block(
?? [a: 0, b: 1, r: isqrt(d), v: [ ], n],
?? if r * r = d then return([r]),
?? do (
????? n: quotient(a + r, b),
????? v: endcons(n, v),
????? a: b * n - a,
????? b: quotient(d - a * a, b),
????? if b = 1 then return(endcons(a + r, v))
?? )
)$
Let's compare results when cf(sqrt(n)) does not yield an integer (here again with cflength=1)
?? a: setify(sublist(1 .. 1000, lambda([n], cf_sqrt(n) # cf(sqrt(n)))))$
?? b: setify(sublist(1 .. 1000, lambda([n], integerp(cf(sqrt(n))))))$
?? setdifference(a, b);
?? {125,216,343,500,864,1000}
For these numbers, the above implementation and Maxima differ. For instance,
?? cf(sqrt(125));
?? [11,4]
?? cf_sqrt(125);
?? [11,5,1,1,5,22]
According to WolframAlpha, the latter is correct, see
http://www.wolframalpha.com/input/?i=continued+fraction+expansion+of+sqrt%28125%29
??
I'm quite puzzled by these results. Is there simply a bug ?
Note: I was just playing with Project Euler problem 64 ;-)
Jean-Claude Arbaut
___________________________________________________________
Retrouvez les scandales de Cannes et les plus belles images du Festival sur Voila.fr http://actu.voila.fr/evenementiel/festival-de-cannes/