41 lines
1.0 KiB
Plaintext
41 lines
1.0 KiB
Plaintext
func tonelli(n, p) {
|
|
legendre(n, p) == 1 || die "not a square (mod p)"
|
|
var q = p-1
|
|
var s = valuation(q, 2)
|
|
s == 1 ? return(powmod(n, (p + 1) >> 2, p)) : (q >>= s)
|
|
var c = powmod(2 ..^ p -> first {|z| legendre(z, p) == -1}, q, p)
|
|
var r = powmod(n, (q + 1) >> 1, p)
|
|
var t = powmod(n, q, p)
|
|
var m = s
|
|
var t2 = 0
|
|
while (!p.divides(t - 1)) {
|
|
t2 = ((t * t) % p)
|
|
var b
|
|
for i in (1 ..^ m) {
|
|
if (p.divides(t2 - 1)) {
|
|
b = powmod(c, 1 << (m - i - 1), p)
|
|
m = i
|
|
break
|
|
}
|
|
t2 = ((t2 * t2) % p)
|
|
}
|
|
|
|
r = ((r * b) % p)
|
|
c = ((b * b) % p)
|
|
t = ((t * c) % p)
|
|
}
|
|
return r
|
|
}
|
|
|
|
var tests = [
|
|
[10, 13], [56, 101], [1030, 10009], [44402, 100049],
|
|
[665820697, 1000000009], [881398088036, 1000000000039],
|
|
[41660815127637347468140745042827704103445750172002, 10**50 + 577],
|
|
]
|
|
|
|
for n,p in tests {
|
|
var r = tonelli(n, p)
|
|
assert((r*r - n) % p == 0)
|
|
say "Roots of #{n} are (#{r}, #{p-r}) mod #{p}"
|
|
}
|