I recall the method you used as being taught in school as a manual way of finding square roots. This was before the era of hand-held calculators. I doubt that it's still being taught in school today. Anyway, you did a remarkable job of converting the method to computer code. I compared your results to those of mine and found no discrepancies other than the exact number of digits displayed, but that's of no consequence.
I used the continued fraction method which was also covered in the Wikipedia article you cited and also in this one:
http://mathworld.wolfram.com/PellEquation.html
The math was a little over my head, but I was able to adapt the general description to code. The independent functions, add2$ and subtract2$ are not part of the solution. They are used to add and subtract the big integers generated by the algorithm.
DIM p(1000) AS SINGLE
DIM q(1000) AS SINGLE
DIM a(1000) AS SINGLE
DIM h(1000) AS STRING
DIM k(1000) AS STRING
INPUT "Number"; tnum
IF SQR(tnum) - INT(SQR(tnum)) = 0 THEN PRINT SQR(tnum), "Perfect Square": END
a(0) = INT(SQR(tnum))
p(1) = a(0)
q(1) = tnum - a(0) ^ 2
a(1) = INT((a(0) + p(1)) / q(1))
h(0) = LTRIM$(STR$(a(0)))
h(1) = LTRIM$(STR$(a(0) * a(1) + 1))
k(0) = "1"
k(1) = LTRIM$(STR$(a(1)))
n = 1
DO
n = n + 1
p(n) = a(n - 1) * q(n - 1) - p(n - 1)
q(n) = (tnum - p(n) ^ 2) / q(n - 1)
a(n) = INT((a(0) + p(n)) / q(n))
FOR i = 1 TO a(n)
h(n) = add2(h(n), h(n - 1))
k(n) = add2(k(n), k(n - 1))
NEXT i
h(n) = add2(h(n), h(n - 2))
k(n) = add2(k(n), k(n - 2))
LOOP UNTIL LEN(k(n)) > 51
temp$ = k(n)
DO
m = 0
DO UNTIL (LEN(h(n)) < LEN(temp$)) OR (LEN(h(n)) = LEN(temp$) AND h(n) < temp$)
h(n) = subtract2(h(n), k(n))
m = m + 1
LOOP
ans$ = ans$ + LTRIM$(STR$(m))
h(n) = h(n) + "0"
LOOP UNTIL LEN(ans$) = 100
dp = LEN(STR$(a(0))) - 1
PRINT LEFT$(ans$, dp) + "." + MID$(ans$, dp + 1)
END
FUNCTION add2$ (a1$, a2$)
alen = LEN(a1$) - LEN(a2$)
IF alen > 0 THEN
a2$ = STRING$(alen, "0") + a2$
ELSE
a1$ = STRING$(ABS(alen), "0") + a1$
END IF
FOR b = LEN(a1$) TO 1 STEP -1
j = VAL(MID$(a1$, b, 1)) + VAL(MID$(a2$, b, 1)) + cry
cry = j \ 10
j$ = LTRIM$(STR$(j MOD 10))
sum$ = j$ + sum$
NEXT b
IF cry THEN sum$ = LTRIM$(STR$(cry)) + sum$
add2$ = sum$
END FUNCTION
FUNCTION subtract2$ (a1$, a2$)
alen = LEN(a1$) - LEN(a2$)
IF alen >= 0 THEN
a2$ = STRING$(alen, "0") + a2$
ELSE
a1$ = STRING$(ABS(alen), "0") + a1$
END IF
FOR b = LEN(a1$) TO 1 STEP -1
d1 = VAL(MID$(a1$, b, 1))
d2 = VAL(MID$(a2$, b, 1))
d = d1 - d2
IF d > 0 THEN j = d1 - d2 - bor: bor = 0
IF d < 0 THEN j = d1 + 10 - d2 - bor: bor = 1
IF d = 0 AND bor = 1 THEN j = d1 + 10 - d2 - bor: bor = 1
IF d = 0 AND bor = 0 THEN j = d1 - d2
j$ = LTRIM$(STR$(j))
dif$ = j$ + dif$
NEXT b
DO
j$ = LEFT$(dif$, 1)
IF j$ = "0" AND LEN(dif$) > 1 THEN dif$ = MID$(dif$, 2)
LOOP UNTIL j$ <> "0" OR LEN(dif$) = 1
subtract2$ = dif$
END FUNCTION
|