BSGS算法_Baby steps giant steps算法（无扩展）详解

Baby Steps-Varsity

Giant Step-Astronauts(May'n・椎名慶治)

阅读时可以听听这两首歌，加深对这个算法的理解。（Baby steps少女时代翻唱过，这个原唱反而不是很有名……Giant Step就比较碉，是一个假面骑士片的插曲，由超碉的May'n和一个人建立的临时组合唱的，怕不怕）

 

这个主要是用来解决这个题：

A^x=B(mod C)(C是质数)，都是整数，已知A、B、C求x。

我在网上看了好多介绍，觉得他们写得都不够碉，我看不懂…于是我也来写一发。

先把x=i*m+j，其中m=ceil(sqrt(C))，（ceil是向上取整）。

这样原式就变为A^(i*m+j)=B(mod C)，

再变为A^j=B*A^(-m*i) (mod C)，

先循环j=0~(C-1)，把（A^j,j)加入hash表中，这个就是Baby Steps（现在可以播放Baby Steps-Varsity了）。

下面我们要做的是枚举等号右边，从hash表中找看看有没有，有的话就得到了一组i j，x=i*m+j，得到的这个就是正确解。

所以，接下来要解决的就是枚举B*A^(-m*i) (mod C)这一步（这就是Giant Step，快放Giant Step-Astronauts(May'n・椎名慶治)）。

A^(-m*i)相当于1/(A^(m*i))，里面有除法，在mod里不能直接用除法，这时候我们就要求逆元。

/*百度百科：

若ax≡1 mod f, 则称a关于模f的乘法逆元为x。也可表示为ax≡1(mod f)。
当a与f互素时，a关于模f的乘法逆元有唯一解。如果不互素，则无解。如果f为素数，则从1到f-1的任意数都与f互素，即在1到f-1之间都恰好有一个关于模f的乘法逆元。
*/
 
然后我们用超碉的exgcd求逆元，exgcd（扩展欧几里德算法）就是在求AB的最大公约数z的同时，求出整数x和y，使xA+yB=z。算法实现就是gcd加几个语句。
然后我们再来看一下exgcd怎么求逆元：

对xA+yB=z，

变成这样xA = z - yB，取B=C（C就是我们要mod的那个）

推导出 xA % C = z %C

只要  z%C==1 时,就可以求出A的逆元x

但用exgcd求完，x可能是负数，还需要这样一下：x=(x%C+C)%C

//--exgcd介绍完毕--

再看我们的题目，

exgcd(A^(m*i) , C)=z，当C是质数的时候z肯定为1，这样exgcd求得的x就是逆元了。

因为x就是A^(m*i)的逆元，P/(A^(m*i))=P*x，所以

B*A^(-m*i) = B/(A^(m*i)) = B*x(mod C)

这样我们的式子A^j=B*A^(-m*i) (mod C)的等号右边就有了，就是B*x，就问你怕不怕！

枚举i，求出右边在hash里找，找到了就返回，无敌！


例题：

http://poj.org/problem?id=2417

代码：

 1 #include<cstdio>
 2 #include<cmath>
 3 #include<iostream>
 4 #include<cstring>
 5 #include<algorithm>
 6 #include<cmath>
 7 #include<map>
 8 #include<set>
 9 using namespace std;
10 #define ll __int64
11 #define usint unsigned int
12 #define RE  freopen("1.in","r",stdin)
13 
14 class hash {
15 public:
16     hash() {
17         memset(a,0xff,sizeof(a));
18     }
19     int locate(ll x) {
20         ll l=x%MOD;
21         while(a[l]!=x&&a[l]!=-1) l=(l+1)%MOD;
22         return l;
23     }
24     void insert(ll x,ll va) {
25         ll l=locate(x);
26         if(a[l]==-1) {
27             a[l]=x;
28             v[l]=va;
29         }
30     }
31     ll get(ll x) {
32         ll l=locate(x);
33         return a[l]==x?v[l]:-1;
34     }
35     void clear() {
36         memset(a,0xff,sizeof(a));
37     }
38 private:
39     static const ll MOD=100007;
40     ll a[MOD+100],v[MOD+100];
41 } S;
42 
43 ll exgcd(ll a,ll b,ll &x,ll &y) {
44     ll t,ret;
45     if (!b) {
46         x=1,y=0;
47         return a;
48     }
49     ret=exgcd(b,a%b,x,y);
50     t=x,x=y,y=t-a/b*y;
51     return ret;
52 }
53 
54 int main() {///A^x=B(modC)    A^j=B*A^(-m*i)(mod C)
55     ll C,A,B;
56     ll m,i,t,D,ans,x,y;
57     while(scanf("%lld%lld%lld",&C,&A,&B)!=EOF) {
58         S.clear();
59         m=ceil(sqrt((double)C));
60         t=1;
61         for(i=0; i<m; i++) { /**One, two, baby steps.Three, four, baby steps.Five, six, baby steps.**/
62             S.insert(t,i);
63             t=t*A%C;
64         }
65         D=1;///此时t=A^m
66         ans=-1;
67         for(i=0; i<m; i++) { /**一歩 Giant Step , 君にとってLittle だとしても**/
68             exgcd(D,C,x,y);///exgcd求逆元，得到x=D^(-i*m)
69             x=((x*B)%C+C)%C;///B*x=B*D^(-i*m)
70             y=S.get(x);
71             //printf("%lld,%lld\n",x,y);
72             if(y!=-1) {
73                 ans=i*m+y;
74                 break;
75             }
76             D=(D*t)%C;///D=t^i,(t=A^m)
77         }
78         if(ans==-1) printf("no solution\n");
79         else printf("%lld\n",ans);
80     }
81     return 0;
82 }
View Code

 

posted @ 2014-07-30 10:24  带鱼Yuiffy  阅读(4606)  评论(1)  编辑  收藏  举报