浅谈同余3(扩展中国剩余定理,扩展BSGS)

发布时间 2023-05-21 13:04:51作者: xuyiyang


距离上一篇已经四个月了,我来填坑了

上一篇:$浅谈同余2(扩展欧几里得,中国剩余定理,BSGS)$

(https://www.cnblogs.com/xyy-yyds/p/17418472.html)

0x50 扩展BSGS $O(\sqrt n)$


【模板】扩展 BSGS/exBSGS

 题目背景

题目来源:SPOJ3105 Mod

题目描述

给定 $a,p,b$,求满足 $a^x≡b \pmod p$ 的最小自然数 $x$ 。

输入格式

每个测试文件中包含若干组测试数据,保证 $\sum \sqrt p\le 5\times 10^6$。

每组数据中,每行包含 $3$ 个正整数 $a,p,b$ 。

当 $a=p=b=0$ 时,表示测试数据读入完全。

输出格式

对于每组数据,输出一行。

如果无解,输出 $No\text{ }Solution$,否则输出最小自然数解。

样例输入 #1

5 58 33
2 4 3
0 0 0

样例输出 #1

9
No Solution

提示

对于 $100\%$ 的数据,$1\le a,p,b≤10^9$ 或 $a=p=b=0$。

2021/5/14 加强 by [SSerxhs](https://www.luogu.com.cn/user/29826)。
2021/7/1 新添加[一组 Hack 数据](https://www.luogu.com.cn/discuss/391666)。

本题$a, p$不互质,所以我们要把$a, p$变得互质
即除去它们的公因数
当然,根据$B\acute{e}zout$定理当$\gcd(a, p) \nmid b$时,方程无解
即求$\frac{a^x}{d} \equiv \frac{b}{d} \pmod {\frac{p}{d}}$ 的答案再加上除的次数就是答案
也可以写成$a ^ {x - cnt} \cdot \frac{a ^ {cnt}}{d} \equiv \frac{b}{d} \pmod {\frac{p}{d}}$
$\frac{a ^ cnt}{d}$要当参数传进去,即代码里的add

注意,本题卡常,此代码要开$C++20 + O2$才能过
我才不会告诉你这道题我卡了一屏TLE呢
$Code:$

#include <iostream>
#include <unordered_map>
#include <cmath>

using namespace std;

typedef long long LL;

inline int read()
{
    register int s = 0, f = 1;
    register char c = getchar();
    for ( ; !isdigit(c); c = getchar())  if (c == '-') f = -1;
    for ( ; isdigit(c); c = getchar())   s = (s << 1) + (s << 3) + (c ^ 48);
    return s * f;
}

inline int qmi(int a, int b, int p)
{
    register int res = 1 % p;
    
    while (b)
    {
        if (b & 1) res = (LL)res * a % p;
        
        b >>= 1;
        a = (LL)a * a % p;
    }
    
    return res;
}

unordered_map<int, int> Hash;
inline int bsgs(int a, int b, int p, int add) //普通BSGS
{
    Hash.clear();
    
    register int t = (int)sqrt(p) + 1;
    register int tmp = 1;
    for (register int j = 0; j < t; j ++ )
    {
        int val = (LL)b * tmp % p;
        Hash[val] = j;
        tmp = (LL)tmp * a % p; //递推减少快速幂的log
    }

    a = qmi(a, t, p);
    tmp = 1;
    for (register int i = 0; i <= t; i ++ )
    {
        int val = (LL)add * tmp % p; //同上
        int j = Hash.find(val) == Hash.end() ? -1 : Hash[val];
        
        if (j >= 0 && i * t - j >= 0) return i * t - j;
        tmp = (LL)tmp * a % p;
    }
    
    return -1;
}

inline int gcd(int a, int b)
{
    return b ? gcd(b, a % b) : a;
}

inline int exbsgs(int a, int b, int p)
{
    a %= p, b %= p;
    if (b == 1 || p == 1) return 0;
    
    register int cnt = 0;
    register int d, add = 1;
    
    while ((d = gcd(a, p)) ^ 1) //消除因子
    {
        if (b % d) return -1; //bezout判无解
        cnt ++, b /= d, p /= d;
        add = (LL)add * a / d % p; //a^cnt / d
        
        if (add == b) return cnt;
    }
    
    int res = bsgs(a, b, p, add);
    if (res == -1) return -1;
    return res + cnt; //别忘记+cnt
}

int main()
{
    int a, b, p, res;
    a = read(), p = read(), b = read();
    while (a || b || p)
    {
        res = exbsgs(a, b, p);
        
        if (res == -1) puts("No Solution");
        else printf("%d\n", res);
        
        a = read(), p = read(), b = read();
    }
    
    return 0;
}

 

0x60 扩展中国剩余定理
求方程组
$$
\begin{cases}
x \equiv {a_1} \pmod {m_1} \\\
x \equiv {a_2} \pmod {m_2} \\\
x \equiv {a_3} \pmod {m_3} \\\
...\\\
x \equiv {a_n} \pmod {m_n} \\\
\end{cases}
$$
的最小整数解,其中$m_1, m_2, ..., m_n$不两两互质

设之前$1 \sim {i - 1}$个方程的解为$x$
前$i - 1$方程的通解可以表示为$x + i * M$, 其中$M = \operatorname{lcm}\\{m_1, m_2, ..., m_{i - 1}\\} $
与$i$个方程合并,可得方程
$x + k * M \equiv a_i \pmod {m_i}$
这就是同余方程,可用$exgcd$求解,然后$x$加上$k * M$就是前$i$个方程的解
要求最小解,所以要对$M$取模

具体可以去看进阶指南和墨染空巨佬的这个(https://www.acwing.com/solution/content/3539/)

以洛谷模板P4777 【模板】扩展中国剩余定理(EXCRT)为例
由于数据范围很大,要炸$long \text{ } long$, 所以要用$int128$
$Code:$

#include <iostream>
#include <cstring>

using namespace std;

typedef __int128 LL;

int n;

LL exgcd(LL a, LL b, LL &x, LL &y) //扩展欧几里得
{
    if (!b)
    {
        x = 1, y = 0;
        return a;
    }
    
    LL d = exgcd(b, a % b, y, x);
    y -= a / b * x;
    
    return d;
}

LL gcd(LL a, LL b)
{
    return b ? gcd(b, a % b) : a;
}

LL lcm(LL a, LL b)
{
    return a / gcd(a, b) * b;
}

int main()
{
    scanf("%d", &n);
    
    LL a1, m1, a2, m2, x, y;
    // x mod m1 = a1
    long long M, A;
    scanf("%lld%lld", &M, &A);
    m1 = M, a1 = A;
    
    for (int i = 2; i <= n; i ++ )
    {
        scanf("%lld%lld", &M, &A);
        m2 = M, a2 = A;
        
        LL d = exgcd(m1, m2, x, y);
        if ((a2 - a1) % d) //bezout判无解
        {
            puts("-1");
            return 0;
        }
        
        x = x * (a2 - a1) / d % (m2 / d);
        if (x < 0) x += m2 / d;
        
        LL mod = lcm(m1, m2);
        a1 = (m1 * x + a1) % mod;
        if (a1 < 0) a1 += mod;
        m1 = mod;
    }
    
    printf("%lld\n", (long long)(a1 % m1));
    
    return 0;
}