03-21-2015, 07:43 AM
The programme CORNACCHIA finds the unique integer solution { x, y } of
x ^ 2 + D * y ^ 2 = p
given D < p & p prime or returns 0 if there is no solution.
The sub-programme SQRTMODP is here
http://www.hpmuseum.org/forum/thread-3448.html
eg For input
{ 23 694283029607 }
the programme returns
{ 829512 16409 } in Ans
& indeed
829512 ^ 2 + 23 * 16409 ^ 2 = 694283029607
Ans►L1:
Ans(1)►Y:
L1(2)►Z:
{-Y,Z}:
RUN SQRTMODP:
IF Ans
THEN
MAX(Z-Ans,Ans)►B:
Z►A:
INT(√Z)►L:
WHILE B>L
REPEAT
A MOD B►R:
B►A:
R►B
END:
√((Z-B^2)/Y)►C:
IF FRAC(Ans)
THEN
0:
ELSE
{B,C}:
END:
ELSE
0:
END:
x ^ 2 + D * y ^ 2 = p
given D < p & p prime or returns 0 if there is no solution.
The sub-programme SQRTMODP is here
http://www.hpmuseum.org/forum/thread-3448.html
eg For input
{ 23 694283029607 }
the programme returns
{ 829512 16409 } in Ans
& indeed
829512 ^ 2 + 23 * 16409 ^ 2 = 694283029607
Ans►L1:
Ans(1)►Y:
L1(2)►Z:
{-Y,Z}:
RUN SQRTMODP:
IF Ans
THEN
MAX(Z-Ans,Ans)►B:
Z►A:
INT(√Z)►L:
WHILE B>L
REPEAT
A MOD B►R:
B►A:
R►B
END:
√((Z-B^2)/Y)►C:
IF FRAC(Ans)
THEN
0:
ELSE
{B,C}:
END:
ELSE
0:
END: