Tonelli-Shanks algorithm
Translation of the Wikipedia pseudocode, heavily influenced by Sidef and Python.
# Legendre operator (πβπ)
sub infix:<β> (Int \π, Int \π where π.is-prime && (πΒ != 2)) {
given π.expmod( (π-1) div 2, π ) {
when 0 { 0 }
when 1 { 1 }
default { -1 }
}
}
sub tonelli-shanks ( \π, \π where (πβπ) > 0 ) {
my $π = π - 1;
my $π = 0;
$π +>= 1 and $π++ while $πΒ %% 2;
return π.expmod((π+1) div 4, π) if $π == 1;
my $π = ((2..π).first: (*βπ) < 0).expmod($π, π);
my $π
= π.expmod( ($π+1) +> 1, π );
my $π‘ = π.expmod( $π, π );
while ($π‘-1)Β % π {
my $b;
my $π‘2 = $π‘Β²Β % π;
for 1 .. $π {
if ($π‘2-1)Β %% π {
$b = $π.expmod(1 +< ($π-1-$_), π);
$π = $_;
last;
}
$π‘2 = $π‘2Β²Β % π;
}
$π
= ($π
* $b)Β % π;
$π = $bΒ²Β % π;
$π‘ = ($π‘ * $π)Β % π;
}
$π
;
}
my @tests = (
(10, 13),
(56, 101),
(1030, 10009),
(1032, 10009),
(44402, 100049),
(665820697, 1000000009),
(881398088036, 1000000000039),
(41660815127637347468140745042827704103445750172002,
100000000000000000000000000000000000000000000000577)
);
for @tests -> ($n, $p) {
try my $t = tonelli-shanks($n, $p);
say "No solution for ({$n}, {$p})." and next ifΒ !$t or ($tΒ² - $n)Β % $p;
say "Roots of $n are ($t, {$p-$t}) mod $p";
}Output:
Last updated
Was this helpful?