离散对数学习:baby_step_giant_step与ex_baby_step_giant_step算法

对于问题 b ^ x == N % mod 的问题，知道b,N,mod，求最小指数 x；

1.对于mod为素数的时候用朴素baby_step_giant_step算法，算法过程如下：

思路： 令x = i * m + j,m为ceil(sqrt(double(mod)));

那么的话等式为 b ^ j == N * b ^ ( - m * i) % mod,要求x，那么就是把满足这个方程的 i , j 求出来，i,j的范围为m = ceil(sqrt(double(mod))),在范围内对左边的打哈希表，这是baby_step操作；然后对右边进行i递增，使得b是以-m的指数进行递增，所以要先求出b%mod的逆元，可以得出(b ^(-m));进行操作的时候如果满足这个等式的话直接说明i,j找到了，然后 x = i * m + j;

复杂度为mlogm

poj2417 

题意：mod为素数，对于问题 b ^ x == N % mod 的问题，知道b,N,mod，求最小指数 x；

思路： mod为素数的话求b的逆元可以用快速幂pow_mod（b,mod - 2,mod）,后者用拓展欧几里得；其他一致

用hash表进行求解是标准的解法：

#include<iostream>#include<cstdio>#include<cstring>#include<cmath>using namespace std;const int maxn = 100000 + 10;#define clr(x,y) memset(x,y,sizeof(x))typedef long long ll;//欧几里得 ll gcd(ll x,ll y){return y ? gcd(y,x % y) : x;}//扩展欧几里得 void ex_gcd(ll a,ll b,ll & d,ll & x,ll & y){if(b == 0){d = a; x = 1;y = 0;return ;}ex_gcd(b,a % b,d,y,x);y -= a / b * x;}//快速幂取模 ll pow_mod(ll x,ll n,ll mod_val){ll ans = 1;ll t = x % mod_val;while(n){if(n & 1){ans = ans * t % mod_val;}n >>= 1;t = t * t % mod_val;}return ans;}//逆元 ll Inv(ll a, ll mod_val){ll x,y,d;ex_gcd(a,mod_val,d,x,y);return d == 1 ? (x % mod_val + mod_val) % mod_val : -1;}//大数相乘取模 ll Mul(ll a,ll b,ll mod_val){a %= mod_val;b %= mod_val;ll ans = 0;ll t = a;while(b){if(b & 1){ans = (ans + t) % mod_val;}b >>= 1;t <<= 1;}return ans;}//哈希表const ll hash_mod = 9876543;ll key[hash_mod],val[hash_mod];int head[hash_mod],next[hash_mod];struct Hash{int tot;void Init(){memset(head,-1,sizeof(head));tot = 0;}void inserts(ll x,ll y){ll k = x % hash_mod;key[tot] = x;val[tot] = y;next[tot] = head[k];head[k] = tot ++;}ll finds(ll x){ll k = x % hash_mod;for(ll i = head[k]; i != -1; i = next[i]){if(key[i] == x){return val[i];}}return -1;}}hs; //求解a^x = b %(mod)已知其他三个，求x; ll baby_giant(ll a,ll b,ll mod_val){hs.Init();ll m = ceil(sqrt((double)mod_val));ll cur = 1;//babyfor(ll j = 0; j < m; j ++){if(hs.finds(cur) == -1){hs.inserts(cur,j);}cur = cur * a % mod_val;}//a^-m %mod_val,这是mod_val为素数的时候的求法 ll _Am = pow_mod(pow_mod(a,mod_val - 2,mod_val),m,mod_val);cur = 1;//giantfor(int i = 0; i < m; i ++){ll j = hs.finds(cur * b % mod_val);if(j != -1){return i *m + j; }cur = cur * _Am % mod_val;}return -1;}//求解a ^x = b %mod ，已知其他三个，求a; ll a,mod,b;int main(){int T = 1; while( ~ scanf("%I64d%I64d%I64d",&mod,&a,&b)){//cout << "case"<<T ++ << ":" <<endl;ll ans = baby_giant(a,b,mod);if(ans != -1)cout << ans << endl;else cout << "no solution" << endl; }return 0;}

下面是用其他的一些方法进行求解的：

map记录位置时间长一些：

#include<iostream>#include<cstdio>#include<cstring>#include<cmath>#include<algorithm>#include<cmath>#include<map>using namespace std;#define clr(x,y) memset(x,y,sizeof x) const int maxn = 1000000 + 10;const double eps = 10e-10;const int mod = 10e6 + 10; typedef long long ll;void ex_gcd(ll a,ll b,ll &d,ll &x,ll &y){if(b == 0){d = a; x = 1;y = 0;return ;}ex_gcd(b,a % b,d,y,x);y -= a/b * x;}ll pow_mod(ll x,ll y,ll mod_val){ll ans = 1;ll t = x % mod_val;while(y){if(y & 1){ans = ans * t % mod_val;}y >>= 1;t = (t * t) % mod_val;}return ans;}ll inv(ll a,ll mod_val){ll x,y,d;ex_gcd(a,mod_val,d,x,y);//x 的通解为c * x + b * k;return d == 1 ? (x % mod_val + mod_val) % mod_val : -1;}ll baby_gaint(ll b,ll n,ll p){ll m = (ll)ceil(sqrt((double)p));ll ans = 1;//baby_stepmap<int,int>ms;//哈希表 for(ll j = 0; j < m; j++){if(!ms[ans])ms[ans] = j + 1;ans = ans * b % p;}//求逆元:当mod为素数的时候，a % mod 的逆元为pow_mod(a,mod - 2,mod)，后者用ex_gcd求; ll invs = pow_mod(b,p - 2,p) // 或者inv(b,p);// b^(-m) %modll _bm = pow_mod(invs,m,p);//gaint_stepans = 1;//大步走m,然后找到与之相等的baby对应的j值 for(ll i = 0; i < m; i ++){if(ms[ans * n % p])//j = ms[ans * n % p]{return i * m + ms[ans * n % p] - 1;}ans = ans * _bm %p;}return -1;}ll b,n,p;int main(){while( ~ scanf("%I64d%I64d%I64d",&p,&b,&n)){ll ans = baby_gaint(b,n,p);if(ans == -1){puts("no solution");continue;}printf("%I64d\n",ans);}return 0;}

用数组记录二分的话实践短很多：

PS:写二分的时候记得找到值不能立即返回，记录下来，可能还有更小的j,

PS： 不知用map建的表会超时，结果用数组记录，二分查找，时间节省了这么多,实践证明不要用map去记录下标，那要很费时间，还是用数组吧；

#include<iostream>#include<cstdio>#include<cstring>#include<cmath>#include<algorithm>#include<cmath>using namespace std;#define clr(x,y) memset(x,y,sizeof x) const int maxn = 1000000 + 10;const double eps = 10e-10;const int mod = 10e6 + 10; typedef long long ll;ll prime[maxn];struct Baby{ll b,j;}baby[maxn];bool cmp(const Baby &x,const Baby & y){return x.b != y.b ? (x.b < y.b) : (x.j < y.j);}ll get_prime(ll x){ll cnt = 0;for(ll i = 2; i<= x; i ++){if(x % i == 0){prime[cnt ++] = i;while(x % i == 0)x /= i;}}if(x > 1)prime[cnt ++] = x;return cnt;}void ex_gcd(ll a,ll b,ll &d,ll &x,ll &y){if(b == 0){d = a; x = 1;y = 0;return ;}ex_gcd(b,a % b,d,y,x);y -= x * a/b;}ll pow_mod(ll x,ll y,ll mod_val){ll ans = 1;ll t = x % mod_val;while(y){if(y & 1){ans = ans * t % mod_val;}y >>= 1;t = (t * t) % mod_val;}return ans;}ll inv(ll a,ll mod_val){ll x,y,d;ex_gcd(a,mod_val,d,x,y);return d == 1 ? (x % mod_val + mod_val) % mod_val : -1;}ll finds(ll x,ll m){ll l = 0,r = m- 1;ll ans = -1;while(r >= l){ll mid = (l + r) / 2;if(baby[mid].b > x){if(baby[mid].b == x)//保证找到的j是最小的 ans = baby[mid].j;r = mid - 1;}else{l = mid + 1;}}return ans;}ll baby_gaint(ll b,ll n,ll p){ll m = (ll)ceil(sqrt((double)(p -1)));//cout << m << endl;ll ans = 1;//baby_stepfor(ll j = 0; j < m; j++){baby[j].j = j;baby[j].b = ans;ans = ans * b % p;}sort(baby,baby+m,cmp);//求逆元:当mod为素数的时候，a % mod 的逆元为pow_mod(a,mod - 2,mod)，后者用ex_gcd求; ll invs = pow_mod(b,p - 2,p);// b^(-m) %modll _bm = pow_mod(invs,m,p);//cout << _bm << endl;//gaint_stepans = 1;//大步走m,然后找到与之相等的baby对应的j值 for(ll i = 0; i < m; i ++){ll j = finds(n * ans % p,m);if(j != -1){//cout <<n * ans %p << " " << i << " " << m << " " << j << endl;return i * m + j;}ans = ans * _bm %p;}return -1;}ll b,n,p;int main(){while( ~ scanf("%I64d%I64d%I64d",&p,&b,&n)){ll ans = baby_gaint(b,n,p);if(ans == -1){cout << "no solution" << endl;continue;}cout << ans << endl;}return 0;}

2.当mod不为素数的时候，用ex_baby_step_gaint_step算法；

算法过程如下：

扩展小步大步攻击：A ^ x = B (mod C)

1 ：  i = 0-> 100 if A^i mod C == B return i;///做这一步是第四步有一个n

2 ：  消因子， 将A^x = B mod C 划为 d * A ^ (x – n) = B (modC)

3 ： m = ceil (sqrt(C))

4 ： 将 A ^ (x - n) 化为  A ^ (I * m + j)(原x 化为了 n + I * m + j)这里与小步大步攻击不同

5 ： for j = 0 ->m   hash（j,A ^ j mod  C）

6 ： for i = 0 ->m   AA = B/(d * A ^ m ^i)

7 ：在hash表中查找AA，若有，取AA对应的j，则答案为I * m + j + n；

这里令x - n = i * m + j,则  A ^j = B * A ^(-i * m) * d ^(-1) % C;求解过程就跟朴素的那个baby_step_gaint_step一样了；

poj3243

hash:

#include<iostream>#include<cstdio>#include<cstring>#include<cmath>using namespace std;const int maxn = 100000 + 10;#define clr(x,y) memset(x,y,sizeof(x))typedef long long ll;//欧几里得ll gcd(ll x,ll y){return y ? gcd(y,x % y) : x;}//扩展欧几里得void ex_gcd(ll a,ll b,ll & d,ll & x,ll & y){if(b == 0){d = a; x = 1;y = 0;return ;}ex_gcd(b,a % b,d,y,x);y -= a / b * x;}//快速幂取模ll pow_mod(ll x,ll n,ll mod_val){ll ans = 1;ll t = x % mod_val;while(n){if(n & 1){ans = ans * t % mod_val;}n >>= 1;t = t * t % mod_val;}return ans;}//逆元ll Inv(ll a, ll mod_val){ll x,y,d;ex_gcd(a,mod_val,d,x,y);return d == 1 ? (x % mod_val + mod_val) % mod_val : -1;}//大数相乘取模ll Mul(ll a,ll b,ll mod_val){a %= mod_val;b %= mod_val;ll ans = 0;ll t = a;while(b){if(b & 1){ans = (ans + t) % mod_val;}b >>= 1;t <<= 1;}return ans;}//哈希表const ll hash_mod = 9876543;ll key[hash_mod],val[hash_mod];int head[hash_mod],nexts[hash_mod];struct hash{int tot;void Init(){memset(head,-1,sizeof(head));tot = 0;}void inserts(ll x,ll y){ll k = x % hash_mod;key[tot] = x;val[tot] = y;nexts[tot] = head[k];head[k] = tot ++;}ll finds(ll x){ll k = x % hash_mod;for(ll i = head[k]; i != -1; i = nexts[i]){if(key[i] == x){return val[i];}}return -1;}}hs;//求解a^x = b %(mod)已知其他三个，求x;ll baby_giant(ll a,ll b,ll mod_val){hs.Init();ll m = ceil(sqrt((double)mod_val));ll cur = 1;//babyfor(ll j = 0; j < m; j ++){if(hs.finds(cur) == -1){hs.inserts(cur,j);}cur = cur * a % mod_val;}//a^-m %mod_val,这是mod_val为素数的时候的求法ll _Am = pow_mod(pow_mod(a,mod_val - 2,mod_val),m,mod_val);cur = 1;//giantfor(int i = 0; i < m; i ++){ll j = hs.finds(cur * b % mod_val);if(j != -1){return i *m + j;}cur = cur * _Am % mod_val;}return -1;}ll ex_baby_giant(ll a,ll b, ll mod_val){hs.Init();for(ll i = 1; i <= 100; i ++){if(pow_mod(a,i,mod_val) == b % mod_val)return i;}ll g,d = 1,n = 0;while((g = gcd(a,mod_val)) != 1){if(b % g)return -1;b /=g;mod_val /= g;n ++;d = d * (a / g) % mod_val;}ll m = ceil(sqrt(double(mod_val)));//babyll cur = 1;for(ll j = 0; j < m; j ++){if(hs.finds(cur) == -1){hs.inserts(cur,j);}cur = cur * a % mod_val;}//d^-1ll _d = Inv(d,mod_val);//d ^ -mll _bm = pow_mod(Inv(a,mod_val),m,mod_val);cur = 1;//giantfor(ll i = 0; i < m; i ++){ll j = hs.finds(cur * _d % mod_val * b % mod_val);if(j != -1){return i * m + j + n;}cur = cur * _bm % mod_val;}return -1;}//求解a ^x = b %mod ，已知其他三个，求a;ll a,mod,b;int main(){int T = 1;while( ~ scanf("%I64d%I64d%I64d",&a,&mod,&b)){if(!a && !b && !mod)break;//cout << "case"<<T ++ << ":" <<endl;//if(b >= mod)//{//cout << "No solution" << endl;//continue;//}ll ans = ex_baby_giant(a,b,mod);if(ans != -1)cout << ans << endl;else cout << "No Solution" << endl;}return 0;}

数组二分：

#include<iostream>#include<cstdio>#include<cstring>#include<cmath>#include<algorithm>#include<cmath>#include<map>using namespace std;#define clr(x,y) memset(x,y,sizeof x) const int maxn = 1000000 + 10;const double eps = 10e-10;const int mod = 10e6 + 10; typedef long long ll;struct Baby{ll b;ll j;}baby[maxn];bool cmp(const Baby &x,const Baby & y){return x.b != y.b ? (x.b < y.b) : (x.j < y.j);}void ex_gcd(ll a,ll b,ll &d,ll &x,ll &y){if(b == 0){d = a; x = 1;y = 0;return ;}ex_gcd(b,a % b,d,y,x);y -= a/b * x;}ll pow_mod(ll x,ll y,ll mod_val){ll ans = 1;ll t = x % mod_val;while(y){if(y & 1){ans = ans * t % mod_val;}y >>= 1;t = (t * t) % mod_val;}return ans;}ll inv(ll a,ll mod_val){ll x,y,d;ex_gcd(a,mod_val,d,x,y);//x 的通解为c * x + b * k;return d == 1 ? (x % mod_val + mod_val) % mod_val : -1;}ll gcd(ll x,ll y){return y ? gcd(y,x%y) : x;}ll finds(ll val,ll m){ll l = 0,r = m - 1;ll ans = -1;while(l <= r){ll mid = (l + r) / 2;if(baby[mid].b > val){if(baby[mid].b == val)//保证返回的j是最小的 {ans =  baby[mid].j;}r = mid - 1;}else {l = mid + 1;}}return ans;}ll ex_baby_gaint(ll A,ll B,ll C){for(int i = 0; i <= 100; i ++){if(pow_mod(A,i,C) == B % C)return i;}ll n= 0,d = 1,g;while((g = gcd(A,C)) != 1){if(B % g)return -1;B /= g;C /= g;n ++;d = d * (A / g) % C;}//cout << A << " " << B << " " << C << " " << n << " " << d << endl;ll m = ceil(sqrt((double)C));ll K = 1;for(ll j = 0; j < m; ++ j){baby[j].b = K;baby[j].j = j;K = K * A % C;}sort(baby,baby + m,cmp);K = 1;ll A_invs = pow_mod(inv(A,C),m,C);//A^-m;ll d_inv = inv(d,C);//d ^ -1if(inv(A,C) == -1 || d_inv == -1){return -1;}for(ll i = 0; i < m; ++ i){ll j = finds(K * d_inv % C * B % C,m);if(j != -1){return i * m + j + n;}K = K * A_invs % C;}return -1;}ll b,n,p;int main(){while( ~ scanf("%I64d%I64d%I64d",&b,&p,&n)){if(!p && !b && !n)break;ll ans = ex_baby_gaint(b,n,p);if(ans == -1){puts("No Solution");continue;}printf("%I64d\n",ans);}return 0;}

hdu2815

注意k >= z的时候进行返回找不到

#include<iostream>#include<cstdio>#include<cstring>#include<cmath>#include<algorithm>#include<cmath>#include<map>using namespace std;#define clr(x,y) memset(x,y,sizeof x) const int maxn = 1000000 + 10;const double eps = 10e-10;const int mod = 10e6 + 10; typedef long long ll;struct Baby{ll b;ll j;}baby[maxn];bool cmp(const Baby & x,const Baby & y){return x.b != y.b ? (x.b < y.b) : (x.j < y.j);}void ex_gcd(ll a,ll b,ll &d,ll &x,ll &y){if(b == 0){d = a; x = 1;y = 0;return ;}ex_gcd(b,a % b,d,y,x);y -= a/b * x;}ll pow_mod(ll x,ll y,ll mod_val){ll ans = 1;ll t = x % mod_val;while(y){if(y & 1){ans = ans * t % mod_val;}y >>= 1;t = (t * t) % mod_val;}return ans;}ll inv(ll a,ll mod_val){ll x,y,d;ex_gcd(a,mod_val,d,x,y);//x 的通解为c * x + b * k;return d == 1 ? (x % mod_val + mod_val) % mod_val : -1;}ll gcd(ll x,ll y){return y ? gcd(y,x%y) : x;}ll finds(ll val,ll m){ll l = 0,r = m - 1;ll ans = -1;while(l <= r){ll mid = (l + r) / 2;if(baby[mid].b >= val){if(baby[mid].b == val){ans = baby[mid].j;}r = mid - 1;}else {l = mid + 1;}}return ans;}ll ex_baby_gaint(ll A,ll B,ll C){for(int i = 0; i <= 100; i ++){if(pow_mod(A,i,C) == B % C)return i;}ll n= 0,d = 1,g;while((g = gcd(A,C)) != 1){if(B % g)return -1;B /= g;C /= g;n ++;d = d * (A / g) % C;}ll m = ceil(sqrt((double)C));ll K = 1;for(ll j = 0; j < m; ++ j){baby[j].b = K;baby[j].j = j;K = K * A % C;}sort(baby,baby + m,cmp);K = 1;ll A_invs = pow_mod(inv(A,C),m,C);//A^-m;ll d_inv = inv(d,C);//d ^ -1for(ll i = 0; i < m; ++ i){ll j = finds(K * d_inv % C * B % C,m);if(j != -1){return i * m + j + n;}K = K * A_invs % C;}return -1;}ll b,n,p;int main(){while( ~ scanf("%I64d%I64d%I64d",&b,&p,&n)){if(n >= p){puts("Orz,I can’t find D!");continue;}ll ans = ex_baby_gaint(b,n,p);if(ans == -1){puts("Orz,I can’t find D!");continue;}printf("%I64d\n",ans);}return 0;}

下面给出baby_step_giant_step和ex_baby_step_giant_step的模板：

ll gcd(ll x,ll y){return y ? gcd(y,x % y) : x;}//扩展欧几里得void ex_gcd(ll a,ll b,ll & d,ll & x,ll & y){if(b == 0){d = a; x = 1;y = 0;return ;}ex_gcd(b,a % b,d,y,x);y -= a / b * x;}//快速幂取模ll pow_mod(ll x,ll n,ll mod_val){ll ans = 1;ll t = x % mod_val;while(n){if(n & 1){ans = ans * t % mod_val;}n >>= 1;t = t * t % mod_val;}return ans;}//逆元ll Inv(ll a, ll mod_val){ll x,y,d;ex_gcd(a,mod_val,d,x,y);return d == 1 ? (x % mod_val + mod_val) % mod_val : -1;}//大数相乘取模ll Mul(ll a,ll b,ll mod_val){a %= mod_val;b %= mod_val;ll ans = 0;ll t = a;while(b){if(b & 1){ans = (ans + t) % mod_val;}b >>= 1;t <<= 1;}return ans;}//哈希表const ll hash_mod = maxn * 10;ll key[hash_mod],val[hash_mod];int head[hash_mod],nexts[hash_mod];struct hash{int tot;void Init(){memset(head,-1,sizeof(head));tot = 0;}void inserts(ll x,ll y){ll k = x % hash_mod;key[tot] = x;val[tot] = y;nexts[tot] = head[k];head[k] = tot ++;}ll finds(ll x){ll k = x % hash_mod;for(ll i = head[k]; i != -1; i = nexts[i]){if(key[i] == x){return val[i];}}return -1;}}hs;//求解a^x = b %(mod)已知其他三个，求x;ll baby_giant(ll a,ll b,ll mod_val){hs.Init();ll m = ceil(sqrt((double)mod_val));ll cur = 1;//babyfor(ll j = 0; j < m; j ++){if(hs.finds(cur) == -1){hs.inserts(cur,j);}cur = cur * a % mod_val;}//a^-m %mod_val,这是mod_val为素数的时候的求法ll _Am = pow_mod(pow_mod(a,mod_val - 2,mod_val),m,mod_val);cur = 1;//giantfor(int i = 0; i < m; i ++){ll j = hs.finds(cur * b % mod_val);if(j != -1){return i *m + j;}cur = cur * _Am % mod_val;}return -1;}ll ex_baby_giant(ll a,ll b, ll mod_val){hs.Init();for(ll i = 1; i <= 100; i ++){if(pow_mod(a,i,mod_val) == b % mod_val)return i;}ll g,d = 1,n = 0;while((g = gcd(a,mod_val)) != 1){if(b % g)return -1;b /=g;mod_val /= g;n ++;d = d * (a / g) % mod_val;}ll m = ceil(sqrt(double(mod_val)));//babyll cur = 1;for(ll j = 0; j < m; j ++){if(hs.finds(cur) == -1){hs.inserts(cur,j);}cur = cur * a % mod_val;}//d^-1ll _d = Inv(d,mod_val);//d ^ -mll _bm = pow_mod(Inv(a,mod_val),m,mod_val);cur = 1;//giantfor(ll i = 0; i < m; i ++){ll j = hs.finds(cur * _d % mod_val * b % mod_val);if(j != -1){return i * m + j + n;}cur = cur * _bm % mod_val;}return -1;}