离散对数&&大步小步算法及扩展

bsgs algorithm

a^x≡b(mod n)

大步小步算法，这个算法有一定的局限性，只有当gcd(a,m)=1时才可以用

原理

此处讨论n为素数的时候。

a^x≡b(mod n)（n为素数）

由费马小定理可知，只需要验证0,1,2...n-1是不是解即可，因为a^n-1 = 1mod(n)

算法过程

1、首先求出a⁰，a¹，a²，...，a^m-1 模上n的值是否为b，存储在e[i]中，求出a^m的逆a^-m

2、下面考虑a^m，a^m+1，...，a^2m-1模上n的值是否为b

此时不用一一检查，如果当中有解，相当于存在e[i]，使得e[i] * a^m = b mod(n)

两边乘上a^-m，e[i] = b * a^-m mod(n)，只需要检查存不存在这样的e[i]即可

3、同理，可以递推检查出a^2m - a^3m-1中解的情况

为了方便，把e[i]存储在map<int, int>x中，x[j]表示满足ei =j 的最小下标(因为可能有多个值相同)

 map<ll, int>x;

 ll log_mod(ll a, ll b, ll n)//n为素数

 {

     //if(b >= n)return -1;

     a %= n;b %= n;//注意题目，如果b >= n是不存在解的

     if(a == )

     {

         if(b == )return ;if(b == 1)return 0;

         else return -;

     }

     ll m = ceil(sqrt(n + 0.5));

     ll v = pow(a, n -  - m, n);

     x.clear();

     x[] = m;

     ll e = ;

     for(int i = ; i < m; i++)//计算a的i次方mod n,并存下来

     {

         e = e * a % n;

         if(!x[e])x[e] = i;

     }

     for(int i = ; i < m; i++)//计算a^(im), a^(im+1),...,a^(im+m-1)

     {

         int num = x[b];

         if(num)return i * m + (num == m ?  : num);

         b = b * v % n;

     }

     return -;

 }

扩展：

当n不为素数的时候，如何求解a^x≡b(mod n) ？

转化成gcd(a, n) = 1即可

如何转化呢？

利用公式：ax≡b(mod n)⇔ax/d ≡ b/d (mod n/d),d=gcd(a,n)（这里是a乘上x，上述求解方程为a^x）

每次用一个a和m和b消去gcd(a, n)，消δ次，每次n /= g, b /= g，迭代更新下一个g，直到g=1

最后a^x变成了a^{x -} ^δ*a'，方程变成：a^{x -} ^δ*a' = b' (mod n') 这里a' b' n'都是消因子之后的值

满足：a' * g = a^δ　　b' * g = b　　n' * g = n g = gcd(a^δ, n)

如果某一步g不整除b'，直接返回-1

如果某一步b' = a'，假设此时消因子次数达到cnt次，那就返回cnt

　　因为消因子次数cnt次等价于

　　原方程：a^x-cnta^cnt = b (mod n)

　　消因子后：a^x-cnta' = b' (mod n')

　　若此时b' = a' 那么等式就可以直接消去a' b'，得到：a^x-cnt= 1 (mod n)，解就是cnt

下面求解：a^{x -} ^δ*a' = b' (mod n')

可以求出a'的逆元，然后乘过去，用上述的大步小步算法求解

下面介绍一个技巧不用求逆元。

先把上述式子写成a^x* tmp = b (mod n) 求出x之后加上δ就是解。

设m = sqrt(n +0.5), x = k * m - q (1 <= k <= m, q <= m)

上式可写成tmp * a ^k*m-q = b (mod n)

等价于　　tmp * a ^k*m = b * a^q(mod n)

可以从1-m枚举k，每次求出tmp * a ^k*m 判断是不是存在b*a^q

^{最开始的大步小步算法也可以这样写}

 #include<bits/stdc++.h>

 using namespace std;

 typedef long long ll;

 ll pow(ll a, ll b, ll m)

 {

     ll ans = ;

     a %= m;

     while(b)

     {

         if(b & )ans = (ans % m) * (a % m) % m;

         b /= ;

         a = (a % m) * (a % m) % m;

     }

     ans %= m;

     return ans;

 }

 ll ext_log_mod(ll a, ll b, ll n)

 {

     if(b >= n)return -;//一些特殊情况的判断

     a %= n;

     if(b == )return ;

     //if(n == 1)return -1;

     ll cnt = ;//记录消因子次数

     ll tmp = ;//存当前a'的值

     for(ll g = __gcd(a, n); g != ; g = __gcd(a, n))

     {

         if(b % g)return -;//不能整除

         b /= g; n /= g; tmp = tmp * a / g % n;

         cnt++;

         if(b == tmp)return cnt;

     }

     ll m = sqrt(n + 0.5);

     ll t = b;

     map<ll, ll>Map;//记录b * a ^ i, i

     Map[b] = ;

     for(int i = ; i <= m; i++)

     {

         b = b * a % n;

         Map[b] = i;

     }

     a = pow(a, m, n);

     for(int k = ; k <= m; k++)//枚举k

     {

         tmp = tmp * a % n;//求出tmp*a^(k*m)

         if(Map.count(tmp))return k * m - Map[tmp] + cnt;

     }

     return -;

 }

 int main()

 {

     ll a, b, p;

     while(scanf("%lld%lld%lld", &a, &p, &b) != EOF)

     {

         ll ans = ext_log_mod(a, b, p);

         if(ans == -)printf("Orz,I can’t find D!\n");

         else printf("%lld\n", ans);

     }

     return ;

 }

巴特西

离散对数&&大步小步算法及扩展

最新文章

热门文章