浅谈二次剩余与Cipolla算法
Preface
数论菜鸡来补一手知识黑洞，二次剩余以前OI时期还真一点没了解过，所以先写个板题先
（虽然当初想着反正到时候有数学巨佬队友带我飞，但多学一点总是好的）
二次剩余又俗称模意义下开根，用于求解\(x^2\equiv n\pmod p\)这样的方程
但注意一般情况下我们只讨论当\(p\)为奇素数时的情况，当\(p\)为任意数的情形下问题则变得异常复杂，因此先略去不谈
而求解二次剩余的主流方法有挺多的，我是选择了容易理解且相对较好实现的Cipolla算法
（以下运算若无特殊说明均指在模\(p\)意义下）
前置知识
勒让德符号
首先我们定义勒让德符号\((\frac np)\)，它的取值为：
- 若\(n=0\)，则\((\frac np)=0\)
- 若\(n\)是模\(p\)意义下的二次剩余，则\((\frac np)=1\)
- 若\(n\)是模\(p\)意义下的非二次剩余，则\((\frac np)=-1\)
这东西在这部分内容中没啥实际用处，但在一些求解二次剩余的数量的问题中大有益处，这里就先不展开了
解的数量
考虑给定\(n,p\)时，若二次剩余存在，则解的数量是多少，首先无解和\(n=0\)的情况比较trivial，因此我们主要讨论\(n\ne 0\)的情形
不妨先假设有任意多组解，设其中两个不相等的解为\(x_1,x_2\)
由于\(x_1^2=n=x_2^2\)，移项后有\((x_1-x_2)(x_1+x_2)=0\)，前面的\(x_1-x_2\)不等于\(0\)，那只能\(x_1+x_2\)等于\(0\)了
因此我们知道若二次剩余有解，则它们一定互为相反数
（注意当\(p\)为奇素数情况下这两个解的奇偶性一定不同，因此不会出现两个解相同的情况）
这也就说明了任意一对相反数都对应一个二次剩余，且这些二次剩余是两两不同的
因此我们有一个推论，那就是二次剩余的数量恰为\(\frac{p-1}{2}\)，而其它的非\(0\)数都是非二次剩余，数量也是\(\frac{p-1}{2}\)
当然也很容易得到勒让德符号的一个性质，\(\sum_{n=0}^{p-1} (\frac np)=0\)
欧拉准则
欧拉准则用于快速判断一个数是否为二次剩余，我们先给出其表达式：
考虑它的证明，这里只给出简易流程，首先当\(n=0\)时结论显然成立
否则我们根据费马小定理，得到\(n^{p-1}=1\)，同时又根据二次探测定理，有\(n^{\frac{p-1}{2}}=\pm 1\)
此时如果\(n\)是二次剩余，则至少存在一个\(x\)，使得\(n^{\frac{p-1}2}= x^{2\times \frac{p-1}2}= x^{p-1}=1\)（具体证明可以用原根，但我不太熟只能口胡下）
若\(n\)是非二次剩余，那么\(n^{\frac{p-1}2}\)就只能等于\(-1\)了
Cipolla算法
大致思想
Cipolla算法是一个基于随机的算法，但由于它的命中率很高，所以用起来也不会有问题
它的大致原理就是先随机在\([0,p-1]\)找一个\(a\)使得\(a^2-n\)是非二次剩余，由于上面讲的性质，因此找到这样的\(a\)的概率大约就是\(\frac{1}{2}\)，因此期望两次就可以找到这样的一个\(a\)
然后就要通过一些神奇的操作来通过这个\(a\)算出满足条件的二次剩余了
扩域
考虑定义模运算下的虚数域，即定义一个\(i\)，满足\(i^2=a^2-n\)（注意由于\(a^2-n\)是非二次剩余，因此这样的\(i\)实际上是不存在的）
这样所有数都可以被表示成\(A+Bi\)的形式，其中\(A,B\)是\([0,p-1]\)的整数
那么我们有接下来的一个重要推论，\((a+i)^{p+1}=n\)
证明的话要先证明两个引理，首先是\(i^p=-i\)，这个比较简单，大致过程为：
然后另一个就是\((A+B)^p=A^p+B^p\)，这个其实在费马小定理的一种证明中就出现过（我们离散数学课上讲的就是这种证明）
考虑利用二项式定理展开左边的式子后，由于\(p\)是质数，因此除了\(C_p^0,C_p^p\)外，其它的组合数\(C_p^i\)分子上\(p!\)中的\(p\)因子都无法被消掉，因此模\(p\)之后都一定为\(0\)，最后剩下来的就是上面的东西了
因此对于上面的推论就很好证明了，稍微化一下式子就有：
所以换句话说，\((a+i)^{p+1}\)就是\(n\)模\(p\)的二次剩余的一个解，那么它的另一个解也很容易得到了
但还有一个小问题，\((a+i)^{p+1}\)是一个扩域之后的数，它的所谓“虚部”一定等于\(0\)吗？
幸运的是确实如此，具体可以用反证法证明，这里由于篇幅所限又菜又懒就略去了
因此我们可以很容易地求解二次剩余问题了
代码实现
模板题见：Luogu P5491 【模板】二次剩余
数论题的板子都是只要理解了原理就非常好写的东西，复杂度的话由于要做快速幂所以是\(O(\log p)\)的，随机找的那部分基本就可以看作常数了
#include<cstdio>
#include<iostream>
#include<random>
#define RI register int
#define CI const int&
using namespace std;
int t,n,mod,w;
struct Complex
{
int x,y;
inline Complex(CI X=0,CI Y=0)
{
x=X; y=Y;
}
friend inline Complex operator * (const Complex& A,const Complex& B)
{
return Complex((1LL*A.x*B.x+1LL*A.y*B.y%mod*w%mod)%mod,(1LL*A.x*B.y%mod+1LL*A.y*B.x%mod)%mod);
}
};
inline int quick_pow(int x,int p=mod-2,int mul=1)
{
for (;p;p>>=1,x=1LL*x*x%mod) if (p&1) mul=1LL*mul*x%mod; return mul;
}
inline Complex quick_pow(Complex x,int p=mod-2,Complex mul=Complex(1,0))
{
for (;p;p>>=1,x=x*x) if (p&1) mul=mul*x; return mul;
}
mt19937_64 rnd(random_device{}());
inline int Cipolla(CI n)
{
if (quick_pow(n,(mod-1)/2)==mod-1) return -1; int a;
for (uniform_int_distribution<int> dist(0,mod-1);;)
{
a=dist(rnd); w=(1LL*a*a%mod-n+mod)%mod;
if (quick_pow(w,(mod-1)/2)==mod-1) break;
}
return quick_pow(Complex(a,1),(mod+1)/2).x;
}
int main()
{
//freopen("CODE.in","r",stdin); freopen("CODE.out","w",stdout);
for (scanf("%d",&t);t;--t)
{
scanf("%d%d",&n,&mod); n%=mod;
if (!n) { puts("0"); continue; }
int x=Cipolla(n); if (!~x) puts("Hola!"); else
{
int y=mod-x; if (x>y) swap(x,y);
printf("%d %d\n",x,y);
}
}
return 0;
}
Postscript
这下数论天坑又补上一个，剩下一个阶与原根还是有点太劲了的说，到时候再慢慢看看吧
