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 }