剩余问题

剩余问题

　　解决形如\(x^k \equiv a \pmod p\)的问题。暴力找似乎代价太高了。我们先从简单的入手。

二次剩余

　　若\(x^2 \equiv a \pmod p\)有解，则称\(a\)是\(p\)的二次剩余。否则称其为非二次剩余。

　　定义\(\text{Legendre}\)符号：

\[\left( \frac{a}{p} \right)=\begin{cases} 1 & a是p的二次剩余 \\ -1 & a不是p的二次剩余 \\ 0 & p \mid a \end{cases} \]

　　（长相奇怪，其实只是记号）

　　仅讨论\(p\)是奇素数的情况。此时有：

\[\left( \frac{a}{p} \right) \equiv a^\frac{p-1}{2} \pmod p \]

　　这个可以快速判断其是否为二次剩余。

　　\(Proof\)：

\[当p \mid a时命题显然成立 \]

\[当p \nmid a时，如果a是p的二次剩余，则有x^2\equiv a \pmod p \]

\[\therefore (x^2)^{\frac{p-1}{2}} \equiv a^{\frac{p-1}{2}} \pmod p \]

\[\therefore x^{p-1} \equiv a^{\frac{p-1}{2}} \pmod p \]

\[由费马小定理易知左边同余1，所以右边同余1 \]

\[如果a不是p的二次剩余，则不存在x使得x^2 \equiv a \pmod p成立 \]

\[\because p是奇质数 \]

\[\therefore 一定能找到x \neq y，使得xy \equiv a \pmod p \]

\[而共有\varphi(p)=p-1个数，其可以两两配对使得上式成立（由逆元的唯一性得到） \]

\[\therefore (p-1)! \equiv a^{\frac{p-1}{2}} \pmod p \]

\[又由威尔逊定理得:(p-1)! \equiv -1 \pmod p \]

\[\therefore a^{\frac{p-1}{2}} \equiv -1 \pmod p \]

　　我们证明了所有的二次剩余得到的结果都是\(1\)，所有的非二次剩余得到的结果都是\(-1\)，整除时为\(0\)，而又不存在第四者，所以命题能成立。

　　然后有\(p\)是奇质数下，二次剩余和非二次剩余的数量一样多。这个很显然，让\(1\)到\(p-1\)的平方在模\(p\)意义下构成集合，而我们有：

\[x^2 \equiv (p-x)^2 \pmod p \]

　　所以我们只要考虑\(1\)到\(\frac{p-1}{2}\)之间的平方即可。

　　若有\(x^2 \equiv y^2 \pmod p\)，\(x \not\equiv y\)，则\(x^2-y^2 \equiv 0 \pmod p \Rightarrow p \mid (x-y)(x+y)\)，由于\(p\)是质数，所以有\(p \mid x-y\)或\(p \mid x+y\)而只能\(x \equiv -y\)，发现两数只能在\(\frac{p-1}{2}\)两边，所以得到\(1\)到\(\frac{p-1}{2}\)之间，两两的平方不同，而这样得到的二次剩余只有\(\frac{p-1}{2}\)个，非二次剩余也只有\(\frac{p-1}{2}\)个。

Cipolla

　　该算法是用来解决二次剩余问题的，思想基于构造。注意下文解决的是\(x^2 \equiv n \pmod p\)，右边的是\(n\)而不是\(a\)了。一般解决\(p\)是质数的情况。

　　如果\(n\)是二次非剩余，显然无解；

　　如果\(n\)是二次剩余，我们随机一个\(a\)，使得\(a^2-n\)为非二次剩余。由非二次剩余和二次剩余数量一样多，所以有\(\frac{1}{2}\)的概率取到，期望下\(2\)次即能得到。定义\(\omega = \sqrt{a^2-n}\)，事实上\(\omega\)在整数域中没有数存在，所以类似于定义虚数单位\(i\)一样，存数也先存成\(a+b\omega\)这种形式，之后会消掉这个东西。

　　我们要\(x^2 \equiv a \pmod p\)，先要知道这几个引理。

　　引理1：\((x+\omega)^p \equiv x^p + \omega^p \pmod p\)。

　　\(Proof\)：

\[由二项式定理得：(x+\omega)^p = \sum_{i=0}^p {p \choose i} x^i\omega^{p-i} \equiv x^p+\omega^p \pmod p \]

\[中间若干项项由于组合数系数被模掉了，这样就\text{Q.E.D}了 \]

　　引理2：\(\omega^p \equiv -\omega\)。

　　\(Proof\)：

\[注意\omega可能不存在于整数域中，所以我们不能直接运用费马小定理 \]

\[\omega^p = \omega^{p-1}\omega = (\omega^2)^\frac{p-1}{2}\omega = (a^2-n)^\frac{p-1}{2}\omega \equiv -\omega \pmod p \]

\[这些我们运用的是上面的定义推导出来的 \]

　　结论：\(x=(a+\omega)^\frac{p+1}{2}\)

　　\(Proof\)：

\[(a+\omega)^{p+1} = (a+\omega)^p(a+\omega) = (a^p+\omega^p)(a+\omega)=(a-\omega)(a+\omega)=a^2-\omega^2=a^2-(a^2-n)=n \equiv x^2 \pmod p \]

\[\Rightarrow (a+\omega)^{p+1} \equiv x^2 \pmod p \]

\[\Rightarrow x \equiv (a+\omega)^\frac{p+1}{2} \pmod p \]

\[\text{Q.E.D} \]

　　那这个式子咋求？直接上复数。根据拉格朗日定理，最后虚部一定为\(0\)。不会拉格朗日定理？其实我也不会那我们反证下吧。

　　\(Proof\)：

\[假设存在x=A+B\omega，B \neq 0，使得x^2 = (A+B\omega)^2 \equiv n \pmod p \]

\[那么一定有A^2 + 2AB\omega + B^2\omega^2 = A^2 + 2AB\omega + B^2(a^2-n) \equiv n \pmod p \]

\[\Rightarrow A^2+B^2(a^2-n)-n \equiv -2AB\omega \pmod p \]

\[\because式子左边的虚部为0，表明式子的右边虚部也为0 \]

\[\therefore必须有:-2AB \equiv 0 \pmod p \]

\[\because B \not\equiv 0 \]

\[\therefore A \equiv 0 \]

\[\therefore B^2\omega^2 \equiv n \pmod p \]

\[\therefore \omega^2 \equiv nB^{-2} \pmod p \]

\[n是二次剩余，B^2是二次剩余，所以B^{-2}也应为二次剩余 \]

\[根据\text{Legendre}符号可以推出二次剩余与二次剩余的积为二次剩余，而\omega^2却不是二次剩余 \]

\[矛盾，所以B=0 \]

　　还有一种情况就是\(p=2\)，此时它是质数但不是奇质数，所以特判掉即可。

　　如果\(x\)有解，其实存在两组解，另一组解为\(p-x\)。

namespace Cipo {
#define rd() (1ll*RAND_MAX*rand() + rand()) // 自定义随机函数（范围够大）
	int n, P, a, t;
	struct comp { int x, y; };
	comp operator * (comp a, comp b) {
		return (comp){(1ll*a.x*b.x+1ll*a.y*b.y%P*t)%P, (1ll*a.x*b.y+1ll*a.y*b.x)%P}; // 类似于复数乘法
	}
	int qpow(int a, int b) {
		int res = 1;
		for (int i = a; b; i = 1ll*i*i%P, b >>= 1)
			if (b & 1) res = 1ll*res*i%P;
		return res;
	}
	comp qpow(comp a, int b) { // 定义复数快速幂
		comp res = (comp){1, 0};
		for (comp i = a; b; i = i*i, b >>= 1)
			if (b & 1) res = res*i;
		return res;
	}
	int sqrt(int n, int P) {
		Cipo::n = n, Cipo::P = P;
		if (P == 2) return n; // P=2特判
		if (!n) return 0; // n=0特判
		if (qpow(n, P-1>>1) != 1) return -1; // 非二次剩余直接排除
		while (a = rd() % P, qpow(t = (1ll*a*a-n+P)%P, P-1>>1) != P-1); // 找出一个a
		return qpow((comp){a, 1}, P+1>>1).x; // 答案
	}
}

　　该算法的复杂度：\(\text{O}(\log p)\)。

　　附期望是\(2\)的证明：

\[\sum_{i=1}^{\infty}\frac{1}{2^i}\times i （这个由期望的定义来：概率×次数）\\ =\sum_{i=0}^{\infty}\frac{1}{2^i} + \frac{1}{2}\sum_{i=0}^{\infty}\frac{1}{2^i}+\dotsb - (1+\frac{1}{2}+\frac{1}{4}+\dotsb) \\ =2+1+\frac{1}{2}+\dotsb-2 \\ =2\times 2-2=2 \]

K次剩余

　　定义同上，只不过把\(2\)换成了\(k\)，但无法套用上面的解法求解。这里仍然假定\(p\)是奇质数。

　　发现构造走不通了！没关系，我们学习了原根！

　　众所周知原根遍历了所有与\(p\)互质的数，所以我们用原根替换它们。

　　比如说\(x^k \equiv a \pmod p\)，设\(g\)为原根，\(a=g^s\)，\(x=g^m\)，则原式变成了\((g^m)^k \equiv g^s \pmod p\)，也就是\(g^{km} \equiv g^s \pmod p\)，\(s\)可以通过\(\text{BSGS}\)解得；\(m\)未知。然后由阶的性质得到\(km \equiv s \pmod{\varphi(p)}\)，也就是\(km \equiv s \pmod{p-1}\)，解\(m\)即可。

　　想要所有解？刚刚的同余方程的解在模\(p-1\)的意义下的所有取值即为所有解。事实上题目一般要求找出任意解，正常解同余方程即可。如果上面的任意一步都出了问题（比如说\(s\)无解，或者同余方程没有解），那就没有解了。

// 码量惊人，博主写了140+行。这里只放最核心部分
// 这里用K次剩余解决二次剩余问题
int solve(int k, int a, int P) { // x^k=a (%p)
	int s = BSGS(g, a, P); // 找出s使得g^s=a
	int d = gcd(k, P-1);
	if (s == -1 || s % d) return -1; // BSGS无解或同余方程无解，说明该问题无解
	return qpow(g, 1ll*(s/d)*inv(k/d, (P-1)/d)%P, P); // 解同余方程km=s(%phi(p))
}

int main() {
	init(); // 初始化质数表

	T = read();
	while (T--) {
		N = read(), P = read(); N %= P;
		if (!N) { printf("0\n"); continue; }
		if (P == 2) { printf("1\n"); continue; }
		getG(P); // 得到模P的原根
		int ans = solve(2, N, P); // 解决K次剩余就写成solve(K, N, P)
		if (ans == -1) printf("Hola!\n");
		else printf("%d %d\n", min(ans, P-ans), max(ans, P-ans));
	}

	return 0;
}