Cipolla's algorithm
# Legendre operator (πβπ)
sub infix:<β> (Int \π, Int \π where π.is-prime && (πΒ != 2)) {
given π.expmod( (π-1) div 2, π ) {
when 0 { 0 }
when 1 { 1 }
default { -1 }
}
}
# a coordinate in a Field of p elements
class Fp {
has Int $.x;
has Int $.y;
}
sub cipolla ( Int \π, Int \π ) {
note "Invalid parameters ({π}, {π})"
and return Nil if (πβπ)Β != 1;
my $Ο2;
my $a = 0;
loop {
last if ($Ο2 = ($aΒ² - π)Β % π)βπ < 0;
$a++;
}
# define a local multiply operator for Field coordinates
multi sub infix:<*> ( Fp $a, Fp $b ){
Fp.new:Β :x(($a.x * $b.x + $a.y * $b.y * $Ο2)Β % π),
Β :y(($a.x * $b.y + $b.x * $a.y) Β % π)
}
my $r = Fp.new:Β :x(1), Β :y(0);
my $s = Fp.new:Β :x($a),Β :y(1);
for (π+1) +> 1, * +> 1 ... 1 {
$r *= $s if $_Β % 2;
$s *= $s;
}
return Nil if $r.y;
$r.x;
}
my @tests = (
(10, 13),
(56, 101),
(8218, 10007),
(8219, 10007),
(331575, 1000003),
(665165880, 1000000007),
(881398088036, 1000000000039),
(34035243914635549601583369544560650254325084643201,
100000000000000000000000000000000000000000000000151)
);
for @tests -> ($n, $p) {
my $r = cipolla($n, $p);
say $rΒ ?? "Roots of $n are ($r, {$p-$r}) mod $p"
Β !! "No solution for ($n, $p)"
}Output:
Last updated
Was this helpful?