【数论】同余（六）：高次同余方程：BSGS算法

default111

已于 2022-03-08 14:31:11 修改

阅读量195

点赞数

分类专栏：算法竞赛文章标签：算法 c++ 数学

于 2022-03-04 23:47:29 首次发布

本文链接：https://blog.csdn.net/default012/article/details/123287579

版权

算法竞赛专栏收录该内容

31 篇文章 3 订阅

订阅专栏

同余问题共7part，我的博客链接：

高次同余方程

只考虑： $a^x\equiv b(\mod m)$

为保持统一，放一下手写的哈希（来自kuangbin的模板）

// (x,y)的HASH
#define MOD 76543
int hs[MOD], head[MOD], nxt[MOD], id[MOD], top;
void insert(int x, int y)
{
    int k = x % MOD;
    hs[top] = x, id[top] = y, nxt[top] = head[k], head[k] = top++;
}
int find(int x)
{
    int k = x % MOD;
    for (int i = head[k]; i != -1; i = nxt[i])
        if (hs[i] == x)
            return id[i];
    return -1;
}

BSGS算法

$\gcd(a,m)=1$ ，令 $x=ks-t,s=\lceil \sqrt m\rceil,0\leq t<s,0<k\leq s+1$ （ $s$ 的取值有严格证明这样复杂度最低），则原式变为 $a^{ks}\equiv ba^{t}(\mod m)$ ，枚举可能的 $t$ ，存下 $ba^t$ 的哈希值，再枚举左边的 $k$ 算 $a^{ks}$ 查表。

// 普通BSGS: a^x==b mod m, gcd(a,m)=1
int BSGS(int a, int b, int m)
{
    if(b == 1 || m == 1) // 注意这里的特判
        return 0;
    int s = sqrt(m) + 1;
    memset(head,-1,sizeof(head));
    top=0;
    for (int i = 0, tmp = b; i < s; i++, tmp = tmp * a % m)
        insert(tmp, i);
    int as = qpow(a, s, m); // 快速幂算a^s
    for (int i = 1, tmp = as; i <= s + 1; i++, tmp = tmp * as % m)
    {
        int qry = find(tmp);
        if (qry == -1)
            continue;
        return i * s - qry;
    }
    return -1;
}

扩展BSGS算法

$\gcd(a,m)=d\neq 1$ ，等价变换为 $a^x+my=b$ 有解的条件是 $d\mid b$ 。两边同除 $d$ ，则变为 $\frac{a}{d}a^{x-1}+y\frac{m}{d}=\frac{b}{d}$ 不断检查 $\gcd(\frac{b}{d},a)$ 并除掉，直到为 $1$ 。
记每次除的数乘积为 $D$ ，则原式变为 $\frac{a^k}{D}a^{x-k}+y\frac{m}{D}=\frac{b}{D}$ 即 $a^{x-k}\equiv \frac{b}{D}\frac{D}{a^k}(\mod \frac{m}{D})$ 用BSGS求解即可。

int exBSGS(int a, int b, int m)
{
    if (b == 1 || m == 1)
        return 0;
    int k = 0, atmp = 1;
    while (1)
    {
        int tmp = __gcd(m, a);
        if (tmp == 1)
            break;
        if(b % tmp)
            return -1;
        b /= tmp, m /= tmp, k++;
        atmp = atmp * (a / tmp) % m;
        if (b == atmp)
            return k;
    }
    return BSGS(a, b * inv(atmp, m) % m, m) + k;
}