Simple Continued Fraction of Square Root using Integer Operations Only

Question

I am trying to find out a way to compute the simple continued fraction of a square root. Simple means that the numerators of the expansion is always one. I have an integer square root function already, so I've already eliminated floating point square roots.

I'm following this generalized algorithm(In section: 6.2.2 An algorithm to find a simple CF for any square-root) . I'm interested in simplifying the repeating steps to only involve only integers (including integer division). In step 3 of the table, it seems the exact square root value is needed in order to compute the next value. Is there a way to eliminate this?

I found python code as well for illustration: https://gist.github.com/rubik/1454917/9c9e552de3f662c79dbb04c9776893a4b2c0a100 Line 24 would be the step that seems non trivial to eliminate floats.

Answer 1

Here is $\sqrt{139}$ which was asked in the stackoverflow qestion..This method is integers only, probably exactly what Fermat used, well suited for hand calculations.

Method described by Prof. Lubin at Continued fraction of $\sqrt{67} - 4$

$$ \sqrt { 139} = 11 + \frac{ \sqrt {139} - 11 }{ 1 } $$ $$ \frac{ 1 }{ \sqrt {139} - 11 } = \frac{ \sqrt {139} + 11 }{18 } = 1 + \frac{ \sqrt {139} - 7 }{18 } $$ $$ \frac{ 18 }{ \sqrt {139} - 7 } = \frac{ \sqrt {139} + 7 }{5 } = 3 + \frac{ \sqrt {139} - 8 }{5 } $$ $$ \frac{ 5 }{ \sqrt {139} - 8 } = \frac{ \sqrt {139} + 8 }{15 } = 1 + \frac{ \sqrt {139} - 7 }{15 } $$ $$ \frac{ 15 }{ \sqrt {139} - 7 } = \frac{ \sqrt {139} + 7 }{6 } = 3 + \frac{ \sqrt {139} - 11 }{6 } $$ $$ \frac{ 6 }{ \sqrt {139} - 11 } = \frac{ \sqrt {139} + 11 }{3 } = 7 + \frac{ \sqrt {139} - 10 }{3 } $$ $$ \frac{ 3 }{ \sqrt {139} - 10 } = \frac{ \sqrt {139} + 10 }{13 } = 1 + \frac{ \sqrt {139} - 3 }{13 } $$ $$ \frac{ 13 }{ \sqrt {139} - 3 } = \frac{ \sqrt {139} + 3 }{10 } = 1 + \frac{ \sqrt {139} - 7 }{10 } $$ $$ \frac{ 10 }{ \sqrt {139} - 7 } = \frac{ \sqrt {139} + 7 }{9 } = 2 + \frac{ \sqrt {139} - 11 }{9 } $$ $$ \frac{ 9 }{ \sqrt {139} - 11 } = \frac{ \sqrt {139} + 11 }{2 } = 11 + \frac{ \sqrt {139} - 11 }{2 } $$ $$ \frac{ 2 }{ \sqrt {139} - 11 } = \frac{ \sqrt {139} + 11 }{9 } = 2 + \frac{ \sqrt {139} - 7 }{9 } $$ $$ \frac{ 9 }{ \sqrt {139} - 7 } = \frac{ \sqrt {139} + 7 }{10 } = 1 + \frac{ \sqrt {139} - 3 }{10 } $$ $$ \frac{ 10 }{ \sqrt {139} - 3 } = \frac{ \sqrt {139} + 3 }{13 } = 1 + \frac{ \sqrt {139} - 10 }{13 } $$ $$ \frac{ 13 }{ \sqrt {139} - 10 } = \frac{ \sqrt {139} + 10 }{3 } = 7 + \frac{ \sqrt {139} - 11 }{3 } $$ $$ \frac{ 3 }{ \sqrt {139} - 11 } = \frac{ \sqrt {139} + 11 }{6 } = 3 + \frac{ \sqrt {139} - 7 }{6 } $$ $$ \frac{ 6 }{ \sqrt {139} - 7 } = \frac{ \sqrt {139} + 7 }{15 } = 1 + \frac{ \sqrt {139} - 8 }{15 } $$ $$ \frac{ 15 }{ \sqrt {139} - 8 } = \frac{ \sqrt {139} + 8 }{5 } = 3 + \frac{ \sqrt {139} - 7 }{5 } $$ $$ \frac{ 5 }{ \sqrt {139} - 7 } = \frac{ \sqrt {139} + 7 }{18 } = 1 + \frac{ \sqrt {139} - 11 }{18 } $$ $$ \frac{ 18 }{ \sqrt {139} - 11 } = \frac{ \sqrt {139} + 11 }{1 } = 22 + \frac{ \sqrt {139} - 11 }{1 } $$

Simple continued fraction tableau:
$$ \begin{array}{cccccccccccccccccccccccccccccccccccccccccc} & & 11 & & 1 & & 3 & & 1 & & 3 & & 7 & & 1 & & 1 & & 2 & & 11 & & 2 & & 1 & & 1 & & 7 & & 3 & & 1 & & 3 & & 1 & & 22 & \\ \\ \frac{ 0 }{ 1 } & \frac{ 1 }{ 0 } & & \frac{ 11 }{ 1 } & & \frac{ 12 }{ 1 } & & \frac{ 47 }{ 4 } & & \frac{ 59 }{ 5 } & & \frac{ 224 }{ 19 } & & \frac{ 1627 }{ 138 } & & \frac{ 1851 }{ 157 } & & \frac{ 3478 }{ 295 } & & \frac{ 8807 }{ 747 } & & \frac{ 100355 }{ 8512 } & & \frac{ 209517 }{ 17771 } & & \frac{ 309872 }{ 26283 } & & \frac{ 519389 }{ 44054 } & & \frac{ 3945595 }{ 334661 } & & \frac{ 12356174 }{ 1048037 } & & \frac{ 16301769 }{ 1382698 } & & \frac{ 61261481 }{ 5196131 } & & \frac{ 77563250 }{ 6578829 } \\ \\ & 1 & & -18 & & 5 & & -15 & & 6 & & -3 & & 13 & & -10 & & 9 & & -2 & & 9 & & -10 & & 13 & & -3 & & 6 & & -15 & & 5 & & -18 & & 1 \end{array} $$

$$ \begin{array}{cccc} \frac{ 1 }{ 0 } & 1^2 - 139 \cdot 0^2 = 1 & \mbox{digit} & 11 \\ \frac{ 11 }{ 1 } & 11^2 - 139 \cdot 1^2 = -18 & \mbox{digit} & 1 \\ \frac{ 12 }{ 1 } & 12^2 - 139 \cdot 1^2 = 5 & \mbox{digit} & 3 \\ \frac{ 47 }{ 4 } & 47^2 - 139 \cdot 4^2 = -15 & \mbox{digit} & 1 \\ \frac{ 59 }{ 5 } & 59^2 - 139 \cdot 5^2 = 6 & \mbox{digit} & 3 \\ \frac{ 224 }{ 19 } & 224^2 - 139 \cdot 19^2 = -3 & \mbox{digit} & 7 \\ \frac{ 1627 }{ 138 } & 1627^2 - 139 \cdot 138^2 = 13 & \mbox{digit} & 1 \\ \frac{ 1851 }{ 157 } & 1851^2 - 139 \cdot 157^2 = -10 & \mbox{digit} & 1 \\ \frac{ 3478 }{ 295 } & 3478^2 - 139 \cdot 295^2 = 9 & \mbox{digit} & 2 \\ \frac{ 8807 }{ 747 } & 8807^2 - 139 \cdot 747^2 = -2 & \mbox{digit} & 11 \\ \frac{ 100355 }{ 8512 } & 100355^2 - 139 \cdot 8512^2 = 9 & \mbox{digit} & 2 \\ \frac{ 209517 }{ 17771 } & 209517^2 - 139 \cdot 17771^2 = -10 & \mbox{digit} & 1 \\ \frac{ 309872 }{ 26283 } & 309872^2 - 139 \cdot 26283^2 = 13 & \mbox{digit} & 1 \\ \frac{ 519389 }{ 44054 } & 519389^2 - 139 \cdot 44054^2 = -3 & \mbox{digit} & 7 \\ \frac{ 3945595 }{ 334661 } & 3945595^2 - 139 \cdot 334661^2 = 6 & \mbox{digit} & 3 \\ \frac{ 12356174 }{ 1048037 } & 12356174^2 - 139 \cdot 1048037^2 = -15 & \mbox{digit} & 1 \\ \frac{ 16301769 }{ 1382698 } & 16301769^2 - 139 \cdot 1382698^2 = 5 & \mbox{digit} & 3 \\ \frac{ 61261481 }{ 5196131 } & 61261481^2 - 139 \cdot 5196131^2 = -18 & \mbox{digit} & 1 \\ \frac{ 77563250 }{ 6578829 } & 77563250^2 - 139 \cdot 6578829^2 = 1 & \mbox{digit} & 22 \\ \end{array} $$

Answer 2

Please see my posting in Stackexchange. Infinite precision. Java BigInteger.

I used a Surd Storage type for infinite precision of the square root of n.

$${\frac{b \sqrt {n} + d}{c} = } $$

$$\frac {b c \sqrt {n} - c d + a_i c^2}{b^2 n - d^2 - (a_i c)^2 + 2 a_i c d} $$

The floor value of the sqrt(n) is only used once. After which the remaining iterations are stored as a surd type. This avoids the rounding errors seen in other algorithms and infinite (memory constrained) resolution can be achieved.

$$a_0 = floor \sqrt{n}$$

$$ a_i = \frac{b_i a_0 + d_i}{c_i} $$

$$ b_{i+1} = b_i c$$

$$ c_{i+1} = b_i^2 n - d_i^2 - (a_i c_i)^2 + 2 a_i c_i d_i $$

$$ d_{i+1} = a_i * (c_i)^2 - c_i * d_i $$

$$ g = gcd(b_{i+1} , c_{i+1} , d_{i+1})$$

$$ b_{i+1}=\frac{b_{i+1}}{g} , c_{i+1}=\frac{c_{i+1}}{g} , d_{i+1}=\frac{d_{i+1}}{g} $$

$$ a_{i+1} = \frac{b_{i+1} a_0 + d_{i+1}}{c_{i+1}} $$

Then for i=0 to i=Maximum_terms, a continued fraction is produced beginning with $$ [a_0;a_1,a_2 ... ,2a_0] $$

I terminate the fraction when the $$a_i = 2 a_0$$ This is the point at which the sequence repeats.

The mathematics was done by Electro World and a very nice video on the mathematics can be found here https://youtu.be/GFJsU9QsytM

A boolean true is returned if a repeat series is found and false if a repeat sequence is not found for the desired precision.

Sample Outputs:

$$ \sqrt{139}=[11;1,3,1,3,7,1,1,2,11,2,1,1,7,3,1,3,1,22] $$ repeat length 18

$$ \sqrt{15}= [3;1,6]$$ repeat length 2

$$ \sqrt{2501}=[50;100]$$ repeat length 1

$$ \sqrt{10807}=[103;1,22,9,2,2,5,4,1,1,1,6,15,1,5,2,1,3,6,34,2,34,6,3,1,2,5,1,15,6,1,1,1,4,5,2,2,9,22,1,206]$$ repeat length 40

A possible two fold speed up would be to look at the palindromic nature of series. In this case 34,2,34. Only half the sequence would need to be determined.

Source Code in Java for BigInteger. Variable precision set by Maximum_terms.

        public static Boolean SquareRootConFrac(BigInteger N) {
BigInteger A,B=BigInteger.ONE,C=B,D=BigInteger.ZERO;
BigInteger A0=N.sqrt(),Bi=B,Ci=C,Di=D,G;
BigInteger TwoA0 = BigInteger.TWO.multiply(A0);
int Frac_Length=0, Maximum_terms=10000; //Precision 10000 terms
String str="";
Boolean Repeat=false, Success=false, Initial_BCD=true;
while(!Repeat) {
    Frac_Length++; Success=!(Frac_Length==Maximum_terms);
    A=((B.multiply(A0)).add(D)).divide(C); Repeat=A.equals(TwoA0)||!Success;`
Bi=B.multiply(C);
Ci=(B.multiply(B).multiply(N)).subtract(D.multiply(D)).subtract(A.multiply(A).multiply(C).multiply(C)).add(BigInteger.TWO.multiply(A).multiply(C).multiply(D));
Di=(A.multiply(C).multiply(C)).subtract(C.multiply(D));
G=Bi.gcd(Ci).gcd(Di);
B=Bi.divide(G);C=Ci.divide(G);D=Di.divide(G);

if(Initial_BCD) {str="["+A+";";System.out.print(str);Initial_BCD=false;}
else            {str=""+A;System.out.print(str);if(!Repeat){str=",";System.out.print(str);}}

}
str="]";System.out.println(str);
str="repeat length ";System.out.print(str);
if(Success) {str=""+(Frac_Length-1);System.out.println(str);}
else        {str="not found";System.out.println(str);}
return Success;
}