退役前挣扎一下 说起来高一寒假就学过一次这玩意,学一次忘一次
其实就是求x^2 \equiv a \pmod p。
有不少算法的样子,比较常用的应该是Cipolla算法,可以求解模奇素数的二次剩余。
首先冷静一下:考虑原根,那么就是要求2x' \equiv a' \pmod \varphi。因为\varphi是偶数,所以一半的a有解,且有两组解。
欧拉判别准则:
考虑a^{\frac{p-1}{2}},如果等于1,说明a存在二次剩余,否则不存在。
证明:用x^2替换a带入,发现等于-1的时候x不存在
于是可以得到这样的算法:
求解x^2 \equiv n \pmod p
找到一个数a,使得a^2 - n非二次剩余,定义i^2 \equiv a^2 - n(类似虚数)。
那么有(a+i)^{\frac{p+1}{2}}及其相反数是n的二次剩余。
证明:
考虑使用二项式定理展开。由Lucas定理,可以发现大部分项的组合数都是0。化简得:
\binom10 \binom10 \sqrt{a^2 - n}^{p+1} + \binom10 \binom11 a \sqrt{a^2 - n}^p + \binom11 \binom10 a^p \sqrt{a^2-n} + \binom11 \binom11 a^{p+1}
大力化简一下,得到这个式子就是n。
于是直接模拟以上算法即可。复杂度期望O(logn)。
例题:SCOI2018 Numazu 的蜜柑
好像没啥好写的啊。。推推式子,发现可以在祖先处解一个方程,并查询子树里有多少符合条件的。启发式合并即可。
这个代码是板子
#include <bits/stdc++.h>
using std::cerr;
using std::endl;
int T, n, P;
inline long long fpow(long long x, int y) {
long long ret = 1;
for ( ; y; y >>= 1, x = x * x % P)
if (y & 1) ret = ret * x % P;
return ret;
}
struct comp { long long x, y; };
long long mul_I;
inline comp operator*(const comp &a, const comp &b) {
return { (a.x * b.x + a.y * b.y % P * mul_I) % P, (a.x * b.y + a.y * b.x) % P };
}
inline comp fpow(comp x, int y) {
comp ret = { 1, 0 };
for ( ; y; y >>= 1, x = x * x)
if (y & 1) ret = ret * x;
return ret;
}
int main() {
srand(time(0)); scanf("%d", &T);
while (T--) {
scanf("%d %d", &n, &P);
if (n == 0) { puts("0"); continue; }
if (fpow(n, (P - 1) / 2) != 1) {
puts("Hola!"); continue;
}
for (long long a = rand() % P; ; a = rand() % P) {
if (!a || fpow((a * a - n + P) % P, (P - 1) / 2) == 1)
continue;
mul_I = (a * a - n + P) % P;
comp qwq = { a, 1 };
qwq = fpow(qwq, (P + 1) / 2);
long long x = qwq.x, y = P - qwq.x;
if (x > y) std::swap(x, y);
printf("%lld %lld\n", x, y);
break;
}
}
return 0;
}