二次剩余
作者:互联网
给出 \(n,p\),求解方程
\[x^2\equiv n~(\operatorname{mod}~p) \]其中保证 \(p\) 是奇素数。
二次剩余数量
我们假设 \(n\) 在模意义下有多个不同的根,其中两个是 \(x_1,x_2\),那么 \(x_1^2\equiv x_2^2\)。
移项加拆括号得:\((x_1-x_2)(x_1+x_2)=0\)。因为 \(x_1\neq x_2\),所以 \(x_1+x_2=0\) 且此方程只有这两个根。
又可以知道,任意一对相反数都对应一个二次剩余,而且这些二次剩余是两两不同的。
也就说二次剩余的数量恰为 \(\frac{p-1}{2}\),其他的非 \(0\) 数都是非二次剩余,数量也是 \(\frac{p-1}{2}\)。
欧拉准则
由费马小定理我们可以知道:\(n^{p-1}\equiv 1~(\operatorname{mod}~p)\)。
由于 \(p\) 是奇素数,所以:\(\left(n^{\frac{p-1}{2}}\right)^2\equiv 1~(\operatorname{mod}~p)\)。
因为 \(1\) 开根后等于 \(1,-1\):\(n^{\frac{p-1}{2}}\equiv \pm 1~(\operatorname{mod}~p)\)。
我们假设 \(n^{\frac{p-1}{2}}\equiv 1~(\operatorname{mod}~p)\),令 \(n=g^k\),\(g\) 为模 \(p\) 的原根。
则有 \(g^{k\frac{p-1}{2}}\equiv 1\),又 \(g^{\varphi(p)}\equiv 1\),所以 \(p-1|k\frac{p-1}{2}\),\(k\) 为偶数。
那么我们可以得出,当 \(n=g^{\frac{k}{2}}\) 时,\(n\) 为二次剩余。
所以 \(n^{\frac{p-1}{2}} \equiv 1\) 与 \(n\) 是二次剩余是等价的,\(n^{\frac{p-1}{2}} \equiv -1\) 与 \(n\) 是非二次剩余是等价的。
Cipolla法求解方程
我们先求出一个 \(a\) 使得 \(a^2-n\) 是非二次剩余。因为它的数量是一半一半,所以随机检验期望 \(2\) 次就能找出这个 \(a\)。
定义 \(i^2=a^2-n\),因为 \(a^2-n\) 是非二次剩余,所以我们求不出来 \(i\) 。
但我们可以类比虚数,定义这一个 \(i\) ,然后所有数变成 \(A+Bi\)。
引理1:\(i^p=-i\)。
证明:\(i^p=i(i^2)^{\frac{p-1}{2}}=i(a^2-n)^{\frac{p-1}{2}}=-i\)
引理2:\((a+b)^p=a^p+b^p\)。
证明:二项式定理展开,消去模 \(p\) 为 \(0\) 的组合数,\((a+b)^p=\sum\limits_{i=0}^p\binom{p}{i}a^ib^{p-i}=a^p+b^p\)
推论:\((a+i)^{p+1}=n\)。
证明:\((a+i)^{p+1}=a^{p+1}+i^{p+1}=a^2-i^2=a^2-(a^2-n)=n\)
因为 \(p+1\) 是偶数,所以 \(x_1=(a+i)^{\frac{p+1}{2}}\) 就是我们求出来的一个解。
至于虚部,可以证明这里虚部为 \(0\),证明自行bdfs。
因为两个解是相反数,那么另一个解即为 \(x_2=p-x_1\)。
至此,我们可以在期望复杂度 \(O(\log p)\) 内求解二次剩余。
代码
我这里写了两个快速幂,其实可以只用写扩域的快速幂。
#include<bits/stdc++.h>
#define pc(x) putchar(x)
#define int long long
using namespace std;
inline int read()
{
int x=0,f=1;char ch=getchar();
while(ch<'0'||ch>'9'){f=ch=='-'?-1:f;ch=getchar();}
while(ch>='0'&&ch<='9'){x=(x<<1)+(x<<3)+(ch^48);ch=getchar();}
return x*f;
}
void write(int x)
{
if(x<0){x=-x;putchar('-');}
if(x>9)write(x/10);
putchar(x%10+48);
}
int T,n,p,i2,x0,x1;
struct CP
{
int x,y;
CP(int xx=0,int yy=0){x=xx;y=yy;}
friend CP operator *(CP x,CP y)
{return CP((x.x*y.x+i2*x.y%p*y.y)%p,(x.x*y.y+x.y*y.x)%p);}
};
CP qpow(CP x,int y)
{
CP res(1,0);
while(y)
{
if(y&1)res=res*x;
x=x*x;y>>=1;
}return res;
}
int qpow(int x,int y)
{
int res=1;
while(y)
{
if(y&1)res=res*x%p;
x=x*x%p;y>>=1;
}return res;
}
bool check(int x){return qpow(x,(p-1)>>1)==1;}
signed main()
{
srand(time(NULL));T=read();
while(T-->0)
{
n=read(),p=read();
if(!n){puts("0");continue;}
if(qpow(n,(p-1)/2)==p-1){puts("Hola!");continue;}
int a=rand()%p;
while(!a||check((a*a-n+p)%p))a=rand()%p;
i2=(a*a%p-n+p)%p;x0=qpow(CP(a,1),(p+1)>>1).x;
x1=p-x0;if(x0>x1)swap(x0,x1);
write(x0),pc(' ');if(x1!=x0)write(x1);
pc('\n');
}return 0;
}
标签:剩余,frac,二次,int,res,CP,equiv 来源: https://www.cnblogs.com/ctldragon/p/16056205.html