二次剩余定义：
在维基百科中，是这样说的：如果q等于一个数的平方模 n，则q为模 n 意义下的二次剩余。例如：x2≡n(mod p)。否则，则q为模n意义下的二次非剩余。
Cipolla算法：一个解决二次剩余强有力的工具，用来求得上式的x的一个算法。
需要学习的数论及数学基础：勒让德符号、欧拉判别准则和复数运算。
勒让德符号：判断n是否为p的二次剩余，p为奇质数。
欧拉定理为xφ(p)≡1(mod p)
当p为素数时，可知φ(p)=p-1，转化为xp-1≡1(mod p)
开根号后为 x(p−1)/2≡±1(mod p)，如果等于1就肯定开的了方，为-1一定开不了。所以x是否为n的二次剩余就用这个欧拉判别准则。
qpow(n,(mod-1)>>1)==mod-1
随机找数a,使得a2−n为复数的虚数单位的平方，即
随机一个数a，然后对a2−n进行开方操作（就是计算他勒让德符号的值），直到他们的勒让德符号为-1为止（就是开不了方为止）。 就是找到一个a满足(a2−n)(p−1)/2=−1。
LL a=1; while(qpow((a*a-n+mod)%mod,(mod-1)>>1)!=mod-1) a=rand()%mod;
建立复数乘法运算（(a+bi)(c+di)=(ac+bd*(-1))+(bc+ad)i）
建立一个类似的域，前面寻找了一个a使(a2−n)(p−1)/2=−1，所以我们定义ω=√(a2−n)。那么现在的ω也像i一样，满足ω2=a2−n=-1
node two(node a,node b)//复数相乘 { node ans; ans.x=(a.x*b.x%mod+a.y*b.y%mod*w%mod)%mod; ans.y=(a.x*b.y%mod+a.y*b.x%mod)%mod; return ans; }
答案=(a+ω)(p+1)/2
根据拉格朗日定理，可以得出虚数处的系数一定为0。
1 node q_pow(node a,LL b){ 2 node res; 3 res.x=1,res.y=0; 4 while(b){ 5 if(b&1)res=two(res,a); 6 a=two(a,a); 7 b>>=1; 8 } 9 return res; 10 }
1 node p; 2 p.x=a,p.y=1,w=(a*a-n+mod)%mod; 3 node ans=q_pow(p,(mod+1)>>1); 4 return ans.x;
2019牛客多校训练营第九场B题为Cipolla算法模板题
1 #include<bits/stdc++.h> 2 #define LL long long 3 using namespace std; 4 const LL mod=1e9+7; 5 struct node 6 { 7 LL x,y; 8 }; 9 LL w; 10 node two(node a,node b)//复数相乘 11 { 12 node ans; 13 ans.x=(a.x*b.x%mod+a.y*b.y%mod*w%mod)%mod; 14 ans.y=(a.x*b.y%mod+a.y*b.x%mod)%mod; 15 return ans; 16 } 17 node q_pow(node a,LL b) 18 { 19 node res; 20 res.x=1,res.y=0; 21 while(b) 22 { 23 if(b&1) 24 res=two(res,a); 25 a=two(a,a); 26 b>>=1; 27 } 28 return res; 29 } 30 LL qpow(LL a,LL b) 31 { 32 LL ans=1; 33 a%=mod; 34 while(b) 35 { 36 if(b&1) 37 ans=ans*a%mod; 38 a=a*a%mod,b>>=1; 39 } 40 return ans; 41 } 42 LL solve(LL n) 43 { 44 if(qpow(n,(mod-1)>>1)==mod-1)//勒让德符号 45 return -1; 46 else if(n==0) 47 return 0; 48 LL a=1;//找随机a 49 while(qpow((a*a-n+mod)%mod,(mod-1)>>1)!=mod-1)//勒让德符号 50 a=rand()%mod; 51 node p; 52 p.x=a,p.y=1,w=(a*a-n+mod)%mod; 53 node ans=q_pow(p,(mod+1)>>1);//求出答案 54 return ans.x; 55 } 56 int main() 57 { 58 int T; 59 scanf("%d",&T); 60 LL q,b,n,x,y,c,t=qpow(2,mod-2); 61 while(T--) 62 { 63 scanf("%lld%lld",&b,&c); 64 q=(b*b-4*c+mod)%mod; 65 n=solve(q); 66 if(n==-1) 67 { 68 printf("-1 -1 "); 69 continue; 70 } 71 x=((b+n)%mod)*t%mod,y=(b-x+mod)%mod; 72 if(x>y) 73 swap(x,y); 74 printf("%lld %lld ",x,y); 75 } 76 return 0; 77 }