53 lines
1.5 KiB
Raku
53 lines
1.5 KiB
Raku
|
# 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";
|
|||
|
|
}
|