【数论】BSGS算法（大步小步算法）

最新推荐文章于 2023-08-29 14:09:28 发布

辞树c

最新推荐文章于 2023-08-29 14:09:28 发布

阅读量8.4k

点赞数 12

分类专栏：数论文章标签： c++ 算法 acm竞赛蓝桥杯

本文链接：https://blog.csdn.net/qq_58207591/article/details/123954286

版权

数论专栏收录该内容

6 篇文章 4 订阅

订阅专栏

引入

来看这样一个数论问题：给定一个质数 p，以及正整数 a，b，求满足同余方程 $a^x \equiv b (mod \; p)$ 的最小非负整数，无满足的则输出 -1。

如果只是简单的枚举，那么要想得出结论，由于循环节最大为 p - 1，就需要枚举 0 ~ p - 1 去验证答案，当 p 的数量级达到 10^9 时，这种枚举显然不能满足算法时间复杂度的需求了。

因此就需要用到接下来要介绍的 BSGS 算法。

BSGS算法

BSGS（baby-step giant-step），即大步小步算法。常用于求解离散对数问题，即上面引入中题目的方程 $a^x \equiv b (mod \; p)$ ，它可以在 $O(\sqrt{p})$ 的时间复杂度下完成求解，不过要求与互质。（不要求与互质的扩展 BSGS 算法留到下一篇再介绍）

来看看 BSGS 是如何做到的：

首先令 $x = A\left \lceil \sqrt{p} \right \rceil - B$ ，其中 $0 \leqslant A, B \leqslant \left \lceil \sqrt{p} \right \rceil$ ，将其代入方程得 $a^{A\left \lceil \sqrt{p} \right \rceil - B} \equiv b (mod\; p)$ ，进行指数分解和移项可以得到 $a^{A\left \lceil \sqrt{p} \right \rceil} \equiv ba^B (mod\; p)$ 。

然后分析一下这个方程，可以发现我们已知参数为和，那么我们就可以通过枚举来计算同余方程右端的所有取值，将它们存入哈希表中（方便后面查找）。

接着对于方程左端，同样可以通过枚举 A 得到所有值，如果在得到的这些值中找到与右端的值相同的 A，就可以通过 $x = A\left \lceil \sqrt{p} \right \rceil - B$ 得到方程的解了。

这样由于 A 和 B 的范围为 $0 \leqslant A, B \leqslant \left \lceil \sqrt{p} \right \rceil$ ，枚举 A 和 B 的时间复杂度就为 $O(\sqrt{p})$ 。

（对于引入中的问题，只要在枚举 A 的过程中发现值与右端中某个值相等，那么就可以直接返回，因为在 A 的增大过程中，是逐渐增大的，因此第一次找到相等情况的就是答案）

BSGS模板

ll bsgs(ll a, ll b, ll mod){
    map<ll, ll> mp;
    ll cur = 1, t = sqrt(mod) + 1;
    for(int B = 1; B <= t; B++){
        cur = cur * a % mod;
        mp[b * cur % mod] = B;
    }
    ll now = cur;
    for(int A = 1; A <= t; A++){
        if(mp[now]) return (ll)A * t - mp[now];
        now = now * cur % mod;
    }
    return -1;
}

例题

1.可爱的质数

题目链接：[TJOI2007] 可爱的质数/【模板】BSGS - 洛谷

解题思路：直接套用上面模板即可。

AC代码：

#include <bits/stdc++.h>
#define lowbit(x) (x & -x)
#define mid (l + r >> 1)
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int, int> PII;
typedef pair<ll, ll> PLL;
const int INF = 0x3f3f3f3f;
const ll LLF = 0x3f3f3f3f3f3f3f3f;
const int mod = 1000000007;

ll bsgs(ll a, ll b, ll mod){
    map<ll, ll> mp;
    ll cur = 1, t = sqrt(mod) + 1;
    for(int B = 1; B <= t; B++){
        cur = cur * a % mod;
        mp[b * cur % mod] = B;
    }
    ll now = cur;
    for(int A = 1; A <= t; A++){
        if(mp[now]) return (ll)A * t - mp[now];
        now = now * cur % mod;
    }
    return -1;
}

int main(){
    //freopen("input.txt", "r", stdin);
    ios::sync_with_stdio(false);
    cin.tie(NULL);
    cout.tie(NULL);
    int p, b, n;
    cin >> p >> b >> n;
    ll ans = bsgs(b, n, p);
    if(ans == -1) cout << "no solution" << endl;
    else cout << ans << endl;
    return 0;
}

2.计算器

题目链接：[SDOI2011]计算器 - 洛谷

解题思路：

对于第 1 项任务， $y^z \; mod \; p$ 可以直接用快速幂进行求解；

对于第 2 项任务，满足 $xy \equiv z (mod \; p)$ 的最小非负整数，可以转换成求满足 $x \equiv z y^{-1}(mod \; p)$ 的最小非负整数，这样只要求出 $zy^{-1} \; mod\; p$ 就可以得到最小的，在这之前需要进行特判，当是的倍数且不是的倍数时是无解的，除了这种情况外就可以使用费马小定理求逆元进行求解。

对于第 3 项任务，满足 $y^x \equiv z(mod \; p)$ 的最小非负整数，可以使用 bsgs 算法进行求解。不过也需要特判：当是的倍数且不是的倍数时是无解的；当是的倍数且也是的倍数时答案是 1。

AC代码：

#include <bits/stdc++.h>
#define lowbit(x) (x & -x)
#define mid (l + r >> 1)
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int, int> PII;
typedef pair<ll, ll> PLL;
const int INF = 0x3f3f3f3f;
const ll LLF = 0x3f3f3f3f3f3f3f3f;
const int mod = 1000000007;

ll quick_pow(ll a, ll b, ll p){
    ll res = 1;
    while(b){
        if(b & 1) res = res * a % p;
        a = a * a % p;
        b >>= 1;
    }
    return res;
}

ll inv(ll x, ll p){
    return quick_pow(x, p - 2, p);
}

ll bsgs(ll a, ll b, ll mod){
    map<ll, ll> mp;
    ll cur = 1, t = sqrt(mod) + 1;
    for(int B = 1; B <= t; B++){
        cur = cur * a % mod;
        mp[b * cur % mod] = B;
    }
    ll now = cur;
    for(int A = 1; A <= t; A++){
        if(mp[now]) return (ll)A * t - mp[now];
        now = now * cur % mod;
    }
    return -1;
}

ll solve1(ll y, ll z, ll p){
    return quick_pow(y, z, p);
}

ll solve2(ll y, ll z, ll p){
    return z * inv(y, p) % p;
}

ll solve3(ll y, ll z, ll p){
    return bsgs(y, z, p);
}

int main(){
    //freopen("input.txt", "r", stdin);
    ios::sync_with_stdio(false);
    cin.tie(NULL);
    cout.tie(NULL);
    int t, k;
    while(cin >> t >> k){
        while(t--){
            ll y, z, p;
            cin >> y >> z >> p;
            if(k == 1){
                cout << solve1(y, z, p) << endl;
            }
            if(k == 2){
                if(!(y % p) && z % p) cout << "Orz, I cannot find x!" << endl;
                else cout << solve2(y, z, p) << endl;
            }
            if(k == 3){
                if(!(y % p)){
                    if(z % p) cout << "Orz, I cannot find x!" << endl;
                    else cout << 1 << endl;
                }
                else{
                    ll ans = solve3(y, z, p);
                    if(ans == -1) cout << "Orz, I cannot find x!" << endl;
                    else cout << ans << endl;
                }
            }
        }
    }
    return 0;
}

3.随机数生成器

题目链接：[SDOI2013] 随机数生成器 - 洛谷

解题思路：要求解这样一个递推式 $x_{i + 1} \equiv ax_i + b (mod \; p)$ 只能一个一个推，但是由于数据比较大，这样一定是超时的，我们尝试将其转换成这样 $x_{i + 1} + \frac{b}{a - 1} \equiv a(x_i + \frac{b}{a - 1}) (mod \; p)$ ，这不就成了等比数列的递推式了嘛！那么就可以得到 $x_n + \frac{b}{a - 1} \equiv a^{n - 1}(x_1 + \frac{b}{a - 1}) (mod \; p)$ ，最后再移个项 $a^{n - 1} \equiv \frac{x_n + \frac{b}{a - 1}}{x_1 + \frac{b}{a - 1}} (mod \; p)$ ，这样就成了 bsgs 的模板样子了，直接套模板开搞！

但这样推下来，会漏掉一些特殊情况，需要进行特判：

x_n = x_1 时，直接输出 1；

a = 0 时，原式变成 $x_{i + 1} \equiv b (mod \; p)$ ，若 b = x_n ，则第二天一定读到，输出 2；否则一定无法读到该页，输出 -1；

a = 1 时，原式变成 $x_{i + 1} \equiv x_i + b (mod \; p)$ ，若 b = 0 ，则变为常数列，由于 x_n = x_1 已经判过了，那这种情况一定无法读到目标页，输出 -1；否则根据等差数列的性质， $\frac{(x_n - x_1)}{b} + 1$ 就是读到目标页的天数。

AC代码：

#include <bits/stdc++.h>
#define lowbit(x) (x & -x)
#define mid (l + r >> 1)
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int, int> PII;
typedef pair<ll, ll> PLL;
const int INF = 0x3f3f3f3f;
const ll LLF = 0x3f3f3f3f3f3f3f3f;
const int mod = 1000000007;

ll quick_pow(ll a, ll b, ll p){
    ll res = 1;
    while(b){
        if(b & 1) res = res * a % p;
        a = a * a % p;
        b >>= 1;
    }
    return res;
}

ll inv(ll x, ll p){
    return quick_pow(x, p - 2, p);
}

ll bsgs(ll a, ll b, ll p){
    map<ll, ll> mp;
    ll cur = 1, t = sqrt(p) + 1;
    for(int B = 1; B <= t; B++){
        cur = cur * a % p;
        mp[b * cur % p] = B;
    }
    ll now = cur;
    for(int A = 1; A <= t; A++){
        if(mp[now]) return (ll)A * t - mp[now] + 1;
        now = now * cur % p;
    }
    return -1;
}

int main(){
    //freopen("input.txt", "r", stdin);
    ios::sync_with_stdio(false);
    cin.tie(NULL);
    cout.tie(NULL);
    int t;
    cin >> t;
    while(t--){
        int p, a, b, x1, xn;
        cin >> p >> a >> b >> x1 >> xn;
        if(x1 == xn){
            cout << 1 << endl;
            continue;
        }
        if(a == 0){
            if(xn == b) cout << 2 << endl;
            else cout << -1 << endl;
            continue;
        }
        if(a == 1){
            if(b == 0) cout << -1 << endl;
            else cout << ((xn - x1) % p + p) % p * inv(b, p) % p + 1 << endl;
            continue;
        }
        ll fz = (xn + b * inv(a - 1, p) % p) % p;
        ll fm = (x1 + b * inv(a - 1, p) % p) % p;
        ll ans = bsgs(a, fz * inv(fm, p) % p, p);
        cout << ans << endl;
    }
    return 0;
}