BSGS算法+逆元 POJ 2417 Discrete Logging
POJ 2417 Discrete Logging
Time Limit: 5000MS | Memory Limit: 65536K | |
Total Submissions: 4860 | Accepted: 2211 |
Description
Given a prime P, 2 <= P < 231, an integer B, 2 <= B < P, and an integer N, 1 <= N < P, compute the discrete logarithm of N, base B, modulo P. That is, find an integer L such that
B^l==N(mod p)
Input
Read several lines of input, each containing P,B,N separated by a space.
Output
For each line print the logarithm on a separate line. If there are several, print the smallest; if there is none, print "no solution".
Sample Input
5 2 1 5 2 2 5 2 3 5 2 4 5 3 1 5 3 2 5 3 3 5 3 4 5 4 1 5 4 2 5 4 3 5 4 4 12345701 2 1111111 1111111121 65537 1111111111
Sample Output
0 1 3 2 0 3 1 2 0 no solution no solution 1 9584351 462803587
1 /*BSGS算法+逆元*/ 2 这个主要是用来解决这个题: 3 4 A^x=B(mod C)(C是质数),都是整数,已知A、B、C求x。 5 6 我在网上看了好多介绍,觉得他们写得都不够碉,我看不懂…于是我也来写一发。 7 8 先把x=i*m+j,其中m=ceil(sqrt(C)),(ceil是向上取整)。 9 10 这样原式就变为A^(i*m+j)=B(mod C), 11 12 再变为A^j=B*A^(-m*i) (mod C), 13 14 先循环j=0~(C-1),把(A^j,j)加入hash表中,这个就是Baby Steps 15 16 下面我们要做的是枚举等号右边,从hash表中找看看有没有,有的话就得到了一组i j,x=i*m+j,得到的这个就是正确解。 17 18 所以,接下来要解决的就是枚举B*A^(-m*i) (mod C)这一步(这就是Giant Step 19 20 A^(-m*i)相当于1/(A^(m*i)),里面有除法,在mod里不能直接用除法,这时候我们就要求逆元。 21 22 /*百度百科: 23 24 若ax≡1 mod f, 则称a关于模f的乘法逆元为x。也可表示为ax≡1(mod f)。 25 当a与f互素时,a关于模f的乘法逆元有唯一解。如果不互素,则无解。如果f为素数,则从1到f-1的任意数都与f互素,即在1到f-1之间都恰好有一个关于模f的乘法逆元。 26 */ 27 28 然后我们用超碉的exgcd求逆元,exgcd(扩展欧几里德算法)就是在求AB的最大公约数z的同时,求出整数x和y,使xA+yB=z。算法实现就是gcd加几个语句。 29 然后我们再来看一下exgcd怎么求逆元: 30 对xA+yB=z, 31 32 变成这样xA = z - yB,取B=C(C就是我们要mod的那个) 33 34 推导出 xA % C = z %C 35 36 只要 z%C==1 时,就可以求出A的逆元x 37 38 但用exgcd求完,x可能是负数,还需要这样一下:x=(x%C+C)%C 39 40 //--exgcd介绍完毕-- 41 42 再看我们的题目, 43 44 exgcd(A^(m*i) , C)=z,当C是质数的时候z肯定为1,这样exgcd求得的x就是逆元了。 45 46 因为x就是A^(m*i)的逆元,P/(A^(m*i))=P*x,所以 47 48 B*A^(-m*i) = B/(A^(m*i)) = B*x(mod C) 49 50 这样我们的式子A^j=B*A^(-m*i) (mod C)的等号右边就有了,就是B*x,就问你怕不怕! 51 52 枚举i,求出右边在hash里找,找到了就返回,无敌! 53 54 /*---------分割线-----------------------*/ 55 #include<cmath> 56 #include<iostream> 57 using namespace std; 58 #include<cstdio> 59 #include<cstring> 60 #define mod 100007 61 #define ll long long 62 struct hash 63 { 64 ll a[mod+10],v[mod+10]; 65 hash(){memset(a,-1,sizeof(a));} 66 int locate(ll x) 67 { 68 ll l=x%mod; 69 while(a[l]!=x&&a[l]!=-1) l=(l+1)%mod; 70 return l; 71 } 72 void insert(ll x,int i) 73 { 74 ll l=locate(x); 75 if(a[l]==-1) 76 { 77 a[l]=x; 78 v[l]=i; 79 } 80 } 81 int get(ll x) 82 { 83 ll l=locate(x); 84 return (a[l]==x)?v[l]:-1; 85 } 86 void clear() 87 { 88 memset(a,-1,sizeof(a)); 89 } 90 }s; 91 void exgcd(ll a,ll b,ll &x,ll &y) 92 { 93 if(b==0) 94 { 95 x=1; 96 y=0; 97 return ; 98 } 99 exgcd(b,a%b,x,y); 100 ll t=x; 101 x=y; 102 y=t-a/b*y; 103 } 104 int main() 105 { 106 ll p,b,n; 107 while(scanf("%I64d%I64d%I64d",&p,&b,&n)==3) 108 { 109 s.clear(); 110 ll m=ceil(sqrt(p)); 111 ll t=1; 112 for(int i=0;i<m;++i) 113 { 114 s.insert(t,i); 115 t=(t*b)%p; 116 } 117 ll d=1,ans=-1; 118 ll x,y; 119 for(int i=0;i<m;++i) 120 { 121 exgcd(d,p,x,y); 122 x=((x*n)%p+p)%p; 123 y=s.get(x); 124 if(y!=-1) 125 { 126 ans=i*m+y; 127 break; 128 } 129 d=(d*t)%p; 130 } 131 if(ans==-1) 132 printf("no solution\n"); 133 else printf("%I64d\n",ans); 134 } 135 136 return 0; 137 }