Calculating the optimal continued fraction of a quadratic irrational

This program finds the optimal continued fraction expansion as far as the end of the first period of a quadratic irrational x = (u+t√d)/v, where d,t,u,v are integers, d >1 and nonsquare, t and v nonzero:
x=a0+ ε1/a1+ ε2/a2+ ···+ εj/aj+ εj+1/aj+1+ ···+ εj+l-1/aj+k-1+ εj+l/···
where the first least period is printed in boldface. We first convert x to (P0+√d)/Q0 where Q0 divides d - P02.
The algorithm is based on the paper Optimal continued fractions by Wieb Bosma, Indag. Math. 49 (1987) 353-379 (see pseudo-code).
We find the first occurrence of equality of complete quotients ξjj+k, j < k. We know that if sj > |Q0|, this gives either a period or half-period. We then look for an inequality ξj+t ≠ ξj+k+t, 1 ≤ t < k. If no such inequality is detected, the OCF period length = l = 2 x NSCF period length = 2k; otherwise the OCF period length = l = NSCF period length = k.
We then progressively work back, checking to see if aj-1=aj+l-1 and εjj+l etc. to get the first least period.

We output the partial numerators and denominators εk and ak, the complete quotients ξk=(Pk+√d)/Qk and the convergents rk/sk, with the first least period highlighted in boldface. Here ak=[ξk] or ak=[ξk]+1 and 0< vk/(2vk+sk) <1/2.
e=1 prints the bk and uk/vk.

Warning. The algorithm can get slowed down by the way I calculate ak=[ξk+vk/(2vk+sk)]. e.g. ξ0=(943+√57)/1681, which has period-length 2730. The author has written a faster BC OCF algorithm based on the NSCF algorithm.

Enter d (1 < d < 1015 and not a square):
Enter t (nonzero):
Enter u:
Enter v:
Enter e: (0 or 1)

Last modified 12th November 2009
Return to main page