Baby-Step-Giant-Step 很酷的算法

                      Baby-Step-Giant-Step

    BSGS算法用于解决形如:      A  ^  x  ≡  B  (  mod  C  ) 的问题。

   学这个算法前需要具备以下知识:快速幂取模、扩展欧几里得、同余知识、哈希表(也可以用map,不过更耗时).. 

  一.     普通的Baby-Step-Giant-Step

    对于普通的BSGS,使用时有个限制条件就是:C只能是一个素数。当C 不是素数的时候,便要用到扩展BSGS,我们先来看普通的…

    做法:首先令 x = a * m + b ;  m = ceil ( sqrt ( C ) ) ;  [ceil: 上取整] 0 <= a < m , 0 <= b < m ;

    Then 原式等价于 ( ( A ^ a ) ^ m ) * A ^ b ≡ B (mod C );

    Then 我们枚举a ,从0 到 m ,对每个a ,算出 A^a ,放到哈希表或map里(用哈希更高效…这里先用map) ,如 h=A^a; mp[h]=a;然后我们用 D 来表示 ( ( A ^ a ) ^ m ) ,用Y来表示 A^B 原式便可化为 : D * Y ≡ B (mod C);

    Then 我们再次枚举 a 从 0 到 m ,可以逐个算出 D 的值,在有了D值得基础上,运用扩展欧几里得可以算出 Y ,现在查找原本记录的表,查看这个Y 值是否存在,如果存在便可返回 a^m + b ;

    如果找完0 到 m ,依然没有发现可行的解,表示无解,退出。

    

 1 #include<stdio.h>
 2 #include<math.h>
 3 #include<map>
 4 using namespace std;
 5 /*
 6   Baby-Step-Giant-Step:   用于解决形如: A^x = B (mod C) 的问题
 7   <><> 普通的解法只能解决C是素数的情况
 8   首先,令 x= a*m + b ; -- m= Ceil(sqrt(c));  --  Ceil 表示上取整
 9   then 原式等价于 ( (A^a )^m )*(A^b)=B (mod c);
10   then 那么我们可以枚举 a: 0-m ,记录到 (a,A^a) ,用H 来表示: (A^a)^m ;
11   那么 原式等价于 H * A^b = B (mod c) ;
12   then 求出 A^m ; 枚举 a: 0-m ,对于每个值判断有没有 A^b 存在,如果有便退出...
13 */
14 // 这个代码不敢保证正确无误,因为没有找到题来A,不过思路应该是没错
15 //了解了普通的BSGS,可以去看下面那份代码,注释得比较详细...
16 #define EPS 0.0000001
17 typedef long long LL;
18 LL pow(LL a,LL b,LL c)    //快速幂取模
19 {
20   LL ans=1;
21   while(b>0)
22   {
23     if(b&1)
24     {
25       ans=ans*a%c;
26     }
27     a=a*a%c;
28     b>>=1;
29   }
30   return ans;
31 }
32 void exGcd(LL a,LL b,LL &d,LL &x,LL &y)
33 {
34   if(b==0)
35   {
36     d=a;
37     x=1;
38     y=0;
39     return ;
40   }
41   exGcd(b,a%b,d,x,y);
42   LL t=x;
43   x=y;
44   y=t-a/b*y;
45 }
46 LL BSGS(LL A,LL B,LL C)
47 {
48   map<LL,int> mp;
49   // LL m=(LL)ceil(sqrt(C));  据英明神武的某朱说Gcc不支持这个ceil
50   LL m=(LL)(sqrt(C)+1-EPS);
51   for(int i=m-1;i>=0;i--)
52   {
53     LL h=pow(A,i,C);
54     mp[h]=i;
55   }
56   LL h=pow(A,m,C);
57   LL d,x,y;
58   for(int i=0;i<m;i++)
59   {
60     LL k=pow(h,i,C);
61     exGcd(k,C,d,x,y);
62     if(B%d!=0) continue;
63     x=((x*B/d%C)+C)%C;
64     if(mp.find(x)!=mp.end()) return mp[x]+i*m;
65   }
66   return -1;
67 }
68 int main()
69 {
70   LL A,B,C;
71   while(scanf("%I64d%I64d%I64d",&A,&B,&C)!=EOF)
72   {
73     printf("%I64d\n",BSGS(A,B,C));
74   }
75   return 0;
76 }
代码

 

  二.     扩展Baby-Step-Giant-Step

       对于C不是素数的情况,不能照搬上面的做法,但也大同小异。

         下面我以hdu的模板题来介绍一下:

      题目链接: http://acm.hdu.edu.cn/showproblem.php?pid=2815

      题意: 给出3个数,A,C,B 表示 A^x ≡ B (mod C ) ,让你求最小的x ; 标准的模板题了,不过题目有些小坑,可以先看下discuss,自己找有点浪费时间,也没多大收获...

     

  1 #include<stdio.h>
  2 #include<stdlib.h>
  3 #include<math.h>
  4 #include<conio.h>
  5 #define Mod 40007       //取模的大小,哈希表的大小...
  6 #define Max 100000     //存放的总数
  7 #define EPS 0.000000001//精度
  8 typedef long long LL;
  9 class Hash             //手写哈希
 10 {
 11   public:
 12     int hs[Mod];       //哈希值  设定的哈希函数为 原值 % Mod ,所以哈希值有可能是 0 ~ Mod-1
 13     int next[Max];     //链表    存放哈希值相等的一条链,他的大小取决于所有原值的数量
 14     LL S[Max];         //存放原值
 15     int H[Max];        //存放所有哈希值
 16     int sn;            //不同原值的数量
 17     int hn;            //不同哈希值的数量
 18     Hash()             //构造函数: 定义Hash类变量时初始化
 19     {
 20       sn=0;
 21       hn=0;
 22       for(int i=0;i<Mod;i++)
 23         hs[i]=0;
 24     }
 25     void clear()       //清空函数
 26     {
 27       sn=0;
 28       for(int i=0;i<hn;i++)
 29         hs[H[i]]=0;
 30       hn=0;
 31     }
 32     void add(LL s)           //加入
 33     {
 34       int ha=abs(s)%Mod;     //计算哈希值
 35       if(hs[ha]==0)          //如果该哈希值还未出现过
 36       {
 37         H[hn++]=ha;          //将该哈希值记录起来,同时哈希值数量加 1
 38       }
 39       sn++;                  //0 表示结尾,所以从1 开始存,原值数量加 1,特别针对 hs数组
 40       S[sn]=s;               //将原值记录起来
 41       next[sn]=hs[ha];       //原本原值记录位置
 42       hs[ha]=sn;             //最新原值记录位置,如果从0 开始存,就无法判断此时是空还是1个值
 43       //比如:5 和 10 有一样的哈希值 ,并且 5 和 10 先后加入 那么有:
 44       //5 加入: next[1] = 0; hs[5] = 1; hs[5] 是哈希值为5 的头,表示第一个原值在1的位置
 45       //10加入: next[2] = 1; hs[5] = 2; 表示第一个哈希值为5的在2,第二个在1,第三个不存在
 46     }
 47     int find(LL s)           //查找
 48     {
 49       int ha=abs(s)%Mod;     //计算哈希值
 50       int k=hs[ha];          //
 51       while(k!=0)
 52       {
 53         if(S[k]==s) return k;//找到
 54         k=next[k];           //下一个节点
 55       }
 56       return 0;              //表示没找到
 57     }
 58 };
 59 int pow(int a,int b,int c)   //快速幂取模
 60 {
 61   LL ans=1;
 62   a%=c;
 63   while(b>0)
 64   {
 65     if(b&1)
 66     {
 67       ans=(LL)ans*a%c;
 68     }
 69     b>>=1;
 70     a=(LL)a*a%c;
 71   }
 72   return ans;
 73 }
 74 void exGcd(int a,int b,int &d,int &x,int &y)     //扩展欧几里得
 75 {
 76   if(b==0)
 77   {
 78     d=a;
 79     x=1;
 80     y=0;
 81     return ;
 82   }
 83   exGcd(b,a%b,d,x,y);
 84   LL t=x;
 85   x=y;
 86   y=t-a/b*y;
 87 }
 88 int Gcd(int a,int b)            //最大公约数
 89 {
 90   return b?Gcd(b,a%b):a;
 91 }
 92 /*
 93   对于 A^x = B (mod C) 这个式子,如果无法保证C是素数,就得使用扩展Baby-Step-Giant-Step
 94   首先: 有 A^x = C * k + B ; 然后,我们一步步地消去 A 和 C 的公因子
 95   令 D 从 1 开始, 原式就可以化为 D*A^x = C*k+B;
 96                               --> D*A/(A,C)*A^(x-1) = C/(A,C)*k +B/(A,C);
 97                               --> ...
 98                               --> D*(A/(A,C))^(co))*A^(x-co) = C/((A,C)^co)*k + B/((A,C)^co)
 99   上面的那个循环知道 A、C互质时结束,没进行一次,D便等于 D*A/(A,C);
100   最后,原式化为 D * A^(x-co) = B (mod C) 如果上面循环中,出现B无法整除的情况就表示无解退出
101   这样,我们便保证了 D 和 A 都和 C 互质,这时候,便可以利用普通的Baby-Step-Giant-Step解决
102   同时,要小心的一点就是,因为这样解出来的最小解是肯定大于co的,但如果真实解小于co,就会出
103   错,所以我们可以事先枚举0 - p 的i,进行特判,这里的p=log(C),因为每次消因子最少也要消去2,
104   所以最多消log(C)次,co最大也就是这个了。
105 */
106 int BSGS(int A,int B,int C)     //扩展 Baby-Step-Giant-Step
107 {
108   Hash hp;
109   for(int i=0;i<50;i++)         //枚举1 - log(C) ,可以稍大一点...
110   {
111     if(pow(A,i,C)==B) return i;
112   }
113   int tmp,co=0;
114   LL D=1,k=1;
115   while((tmp=Gcd(A,C))!=1)      //消因子循环
116   {
117     co++;
118     if(B%tmp!=0) return -1;
119     C/=tmp;
120     B/=tmp;
121     D=D*A/tmp%C;
122   }
123   hp.clear();
124   int m=(int)(sqrt((double)C)-EPS);   //m=ceil(sqrt(C)),这里的ceil表示上取整
125   for(int i=0;i<m;i++,k=k*A%C)
126     if(hp.find(k)==0) hp.add(k);
127   int d,x,y;
128   for(int i=0;i<m;i++)
129   {
130     exGcd(D,C,d,x,y);
131     x=((LL)x*B%C+C)%C;                //不用进行B是否能整除d的判断,因为D和C永远互质
132     int w=hp.find(x);
133     if(w!=0) return w-1+co+i*m;
134     D=D*k%C;                          //D最开始的时候和C互质,同时k也和C互质
135   }
136   return -1;                          //找不到
137 }
138 int main()
139 {
140   int A,B,C;
141   while(scanf("%d%d%d",&A,&C,&B)!=EOF)
142   {
143     if(B>=C)
144     {
145       printf("Orz,I can’t find D!\n");
146       continue;
147     }
148     int t=BSGS(A,B,C);
149     if(t==-1)
150       printf("Orz,I can’t find D!\n");
151     else printf("%d\n",t);
152   }
153   return 0;
154 }

 

posted @ 2015-08-29 15:10  hchlqlz  阅读(1282)  评论(0编辑  收藏  举报