POJ1808 平方（二次）同余方程

POJ1808  给定一个方程 x*x==a(mod p) 其中p为质数 判断是否有解

程序中 MOD_sqr（）返回整数解 无解返回-1

数学证明略 

#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<cmath>
#include<vector>
#include<set>
#include<map>
#include<algorithm>
using namespace std;
typedef long long int LL;

 vector<LL>anss;
 LL ans[1000000];
LL pow(LL n,LL k,LL p);
inline LL findv(LL x,LL p);
LL GCD(LL a,LL b)//Euclid 
{
 return b==0?a:GCD(b,a%b);
}
void extend_GCD(LL a,LL b,LL &x,LL &y)//extend Euclid
{
 if(b==0){x=1;y=0;}
    else {extend_GCD(b,a%b,y,x);y-=a/b*x;} 
}
LL line_mod_equation(LL a,LL b,LL n)//二元一次线性不定方程求解 
{
 LL x,y;
 extend_GCD(a,n,x,y);
 LL d=GCD(a,n);
 if(b%d==0)
 	{
 	 while(x<0)x+=n;
 	 while(x>n)x-=n;
 	 LL a1=x*(b/d)%(n/d);
 	 while(a1-n/d>0)a1-=n/d;
 	 return a1;//此处只返回最小正整数解 若要解序列 在主函数中生成 
	}
 return -1;
}
LL CRT(int a[],int m[],int n)//中国剩余定理 
{
 int M=1;
 for(LL i=0;i<n;i++)
    M*=m[i];
 LL ret=0;
 for(LL i=0;i<n;i++)
	{
	 LL x,y;
	 LL tm=M/m[i];
	 extend_GCD(tm,m[i],x,y);
	 ret=(ret+tm*x*a[i])%M;
	}
 while(ret<0)ret+=M;
 while(ret>M)ret-=M;
 return ret;
}
inline LL findRoot(LL p)
{
 LL s=ceil(sqrt(p-1));
 for(LL i=1;i<p;i++)
 	{
 	 if(pow(i,p-1,p)==1)
 	 	{
 	 	 bool sign=1;
 	 	 for(LL j=2;j<=s;j++)
 	 	 	if(((p-1)%j==0)&&(pow(i,(p-1)/j,p)==1))
 	 	 	  {sign=0;break;}
 	 	 if(sign)return i;
 	 	 	
		}
	}
 return -1;
}
LL pow(LL n,LL k,LL p)
{
 if(k==0)return 1;
 LL w=1;
 if(k&1)w=n%p;
 LL ans=pow(n*n%p,k>>1,p);
 return ans*w%p;
}
inline LL ind(LL x,LL p,LL g)//g^y=x(mod p) return y
{
 static map<LL,LL>tmp;
 LL s=ceil(sqrt(p));
 LL w=pow(g,s,p);
 for(LL r=0;r<=s;++r)tmp.insert(make_pair(pow(g, r, p), r));
 for(LL t=0;t<=s;++t)
 	{
 	 LL gt=pow(w,t,p);
 	 LL anti=findv(gt,p);
 	 if(tmp.count(x*anti%p))return (tmp[x*anti%p]+t*s);
 	 
	}
 return -1;
}
inline LL findv(LL x,LL p)
{
 LL d,t;
 extend_GCD(x,p,d,t);
 while(d<0)d+=p;
 return d;
}
LL MOD_sqr(LL a,LL p)//求解 x^2=a(modp)原理有待证明 
{
 if(p==2)return a%p;
 if(pow(a,(p-1)/2,p)!=1)return -1;
 LL x;
 if(p%4==3)x=pow(a,(p+1)/4,p);
    else
    	{
    	 LL b;
    	 for(b=1;pow(b,(p-1)/2,p)==1;b++);
    	 LL i=(p-1)/2,k=0;
    	 do
    	  {
    	   i>>=1;k>>=1;
    	   if((pow(a,i,p)*pow(b,k,p)+1)%p==0)
    	     {
    	      k+=(p+1)/2;
			 }
			
		  }while(i%2==0);
		 x=pow(a,(i+1)/2,p)*pow(b,k/2,p)%p;
		}
 if(x*2>p)x=p-x;
 return x;
} 
int main()
{
 freopen("t.txt","r",stdin);
 LL np;
 scanf("%lld",&np);
 for(LL ii=1;ii<=np;ii++)
   {
    LL a,p;
    scanf("%lld%lld",&a,&p);
    a=(a%p+p)%p;
    LL ans=MOD_sqr(a,p);
    if(ans==-1)printf("Scenario #%lld:
-1

",ii);
       else printf("Scenario #%lld:
1

",ii);
   }
 return 0;
}