数论：从提高组到提高组

Bc2_Ch1ckenPr1nce · 2026-03-04 20:22:45 · 算法·理论

说明

最近可爱的 MGJ 连上了七天的数论课，然后她写了一篇提高组数论的合集。

由于篇幅原因，大部分例题不放代码。

此文章部分参考他人文章或借他人文章进行优化，链接如下：

1. https://www.luogu.com.cn/problem/solution/P5656；
2. https://www.luogu.com.cn/article/fcv6pwn2；
3. https://www.luogu.com.cn/article/n6if1g3y。

更新日志：

• 2026.3.4 22:35：更新扩展欧几里得（Exgcd）；
• 2026.3.5：更新乘法逆元；
• 2026.3.7：更新欧拉函数；
• 2026.3.8-3.9：更新中国剩余定理（CRT）；
• 2026.3.9-3.10：更新欧拉定理；
• 2026.3.11：更新小步大步算法（BSGS）；
• 2026.3.11 21:36：完善所有内容，整改错别字。

本文耗时 7 天在每天课余时间抽空写完，作者写文章不易，如有不当之处敬请谅解（可以私信作者）。

扩展欧几里得（Exgcd）

Part 0：前置知识

• 基础数学（四则运算，符号定义等）；
• 欧几里得定理：\gcd(a, b) = \gcd(b,\ a \bmod b)；
• 裴蜀定理。

Part 1：扩展欧几里得的定义

扩展欧几里得是一个算法，一般解决如下核心问题：

给定两个整数 a, b，找到整数 x, y 满足：

这个方程叫做裴蜀方程（Bézout's identity）。

它在信息学竞赛中有着广泛应用，如求解一次不定方程或同余方程。

Part 2：裴蜀方程的求解

我们使用递归推导。

我们对 ax + by = \gcd(a, b) 用一次欧几里得定理，得：

将 a \bmod b 转化为 a - b \cdot \left\lfloor \frac{a}{b} \right \rfloor，得：

与原式比较，得到一组特解：

于是，我们找到了上层 x, y 和本层 x_1, y_1 的关系。

最终当 b = 0 时，\gcd(a, 0) = a，方程变为 ax + 0y = a，显然得到一组解：x = 1，y = 0。

::::warning[注意] 信息学大忌：在 b = 0 时耍帅把 y 随便写（如 114514），因为每一层递归时 x, y 都会变化。当递归次数足够多时，x, y 将会超出变量范围或产生某些其他问题。 ::::

下面给出两种模板代码，第一种好理解但代码较长（也长不了多少），新手建议使用第一种：

::::success[Code 1]

ll Exgcd(ll a, ll b, ll &x, ll &y) {
  if (!b) {
    x = 1, y = 0;
    return a;
  }
  ll d = Exgcd(b, a % b, x, y), tmp = x;
  x = y, y = tmp - (a / b) * y;
  return d;
}

::::

::::success[Code 2]

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;
}

::::

Part 3：扩展欧几里得的应用

3.1：解不定方程

常见问题：求不定方程 ax + by = c 的一组解。

首先，先使用裴蜀定理判断是否有解：由裴蜀定理得 \gcd(a, b) \mid (ax + by)，故 \gcd(a, b) \nmid c 时不定方程无解。

然后，用 Exgcd 求出 ax + by = \gcd(a, b) 的一组整数特解 x_0, y_0。那么有：

两边同乘 \cfrac{c}{\gcd(a, b)}，得：

与原式比较，可以得到一组特解：

有了特解，那么就可以推导通解了。

显然有：

令 \Delta x = x - x_0，\Delta y = y - y_0。

则：

展开，得：

又 ax_0 + by_0 = c，则：

通过各种方法，解得：

于是，我们得到通解形式：

其中 k \in \mathbb{Z}。

我们还可以推导出，如果要求最小非负整数解，可以提前将系数 a, b 都除以 \gcd(a, b)。

3.2：解同余方程

常见问题：解形如 ax \equiv c \pmod b 的同余方程。

考虑将此方程尽可能写成 ax + by = c 的形式：

于是就变成了不定方程的形式。

Part 4：例题

P1082 [NOIP 2012 提高组] 同余方程 / P2613【模板】有理数取余

两题均为套模板题。

讲解几个需要注意的点：

• 在 P1082 中，为了避免 x_0 为负数，我们需要先 x \to x \bmod b，再将 x \to x + b。但这可能并不是最小整数解，所以我们需要调整为 x' = (x_0 \bmod b + b) \bmod b；
• 在 P2613 中，非负整数 a, b 的范围达到 10 ^{10001}，所以需要结合快读一边读入一边取模；
• 在 P2613 中，需要判断无解情况。

P1516 青蛙的约会

如果他们相遇，他们初始的位置坐标之差和跳的距离应该在模 l 下同余。

令 k 为跳的次数，则 (n - m)k \equiv x - y \pmod l。

套模板即可，注意 n - m 和 x - y 可能为负数和无解情况。

P5656 【模板】二元一次不定方程 (exgcd)

模板简单，模板题不一定简单。因为此题细节很多，需要一一处理。

下面所有变量源于 Part 3 的应用 1。

那么我们先求最小值。首先，肯定的是用 Exgcd 求出一组特解 x_0, y_0，那么显然 x, y 最小正整数解为 x_{\min} = (x_0 \bmod \Delta x + \Delta x) \bmod \Delta x，y_{\min} = (y_0 \bmod \Delta y + \Delta y) \bmod \Delta y（读者可以简单证明）。

然后就可以求最大值了。最大值怎么求？很显然，对于 x_{\min}, y_{\min}，有 ax_{\min} + by_{\max} = c 和 ax_{\max} + by_{\min} = c。那么 y_{\max} = \cfrac{c - ax_{\min}}{b}，x_{\max} = \cfrac{c - by_{\min}}{a}。

你以为这就结束了？我们还没有处理没有正整数解的情况。当我们限制 x, y > 0 时，可以推导出 k 的取值范围，当 k 的下界大于上界时就没有正整数解。

::::success[Code]

#include <bits/stdc++.h>

using namespace std;

using ll = long long;
using ull = unsigned long long;

const int kMaxN = 2e5 + 10;

ll a, b, c, d, x, y, dx, dy, m, n, xmn, xmx, ymn, ymx;

ll Exgcd(ll a, ll b, ll &x, ll &y) {
  if (!b) {
    x = 1, y = 0;
    return a;
  }
  ll res = Exgcd(b, a % b, y, x);
  y -= a / b * x;
  return res;
}

int main() {
  ios::sync_with_stdio(0), cin.tie(0), cout.tie(0);
  int t;
  for (cin >> t; t; -- t) {
    cin >> a >> b >> c;
    d = Exgcd(a, b, x, y);
    if (c % d) {
      cout << "-1\n";
      continue;
    }
    dx = b / d, dy = a / d;
    x *= c / d, y *= c / d;
    m = ceil(1.0 * (1 - x) / dx), n = floor(1.0 * (y - 1) / dy);
    if (m > n) {
      cout << x + m * dx << ' ' << y - n * dy << '\n';
    } else {
      xmn = (x % dx + dx) % dx;
      if (!xmn) {
        xmn = dx;
      }
      ymn = (y % dy + dy) % dy;
      if (!ymn) {
        ymn = dy;
      }
      xmx = (c - ymn * b) / a, ymx = (c - xmn * a) / b;
      cout << n - m + 1 << ' ' << xmn << ' ' << ymn << ' ' << xmx << ' ' << ymx << '\n';
    }
  }
  return 0; 
}

::::

乘法逆元

Part 0：前置知识

• 基础数学（四则运算，符号定义等）；
• 费马小定理。

Part 1：乘法逆元的定义

1.1：为什么需要乘法逆元？

我们先来探讨一下，为什么需要乘法逆元？

因为，对于加减乘运算，我们都可以通过四则运算和模运算进行取模。但是，对于除法：(a \div b) \bmod m \neq ((a \bmod m) \div (b \bmod m)) \bmod m。

1.2：乘法逆元的定义

对于整数 a 和模数 m（m > 1），如果存在整数 x 满足：

ax \equiv 1 \pmod p

则称 x 是 a 在模 m 意义下的乘法逆元，记作 a^{-1} \equiv x \pmod p。

::::warning[注意] 这里的 a^{-1} 不直接等同于 \frac{1}{a}，必须要在模运算下才可以完全等同。 ::::

1.3：逆元存在的条件

定理：a 在模 p 意义下有乘法逆元的充分必要条件是 \gcd(a, p) = 1。

证明：

• 必要性：如果 \gcd(a, p) > 1，假设存在逆元 x，那么 ax = 1 + kp。由于 \gcd(a, p) \mid a 且 \gcd(a, p) \mid p，所以 \gcd(a, b) \mid 1，矛盾。
• 充分性：如果 \gcd(a, p) = 1，根据 Exgcd，存在 x, y \in \mathbb{Z} 使得 ax + py = 1，于是 ax \equiv 1 \pmod p。

Part 2：乘法逆元的求法

2.1：Exgcd 求逆元

• 原理：将同余方程转化为不定方程；
• 优点：有判断无解的步骤，当 p 不是质数时也可使用；
• 时间复杂度：单个 \mathcal{O}(\log x)。

我们将同余方程 a \cdot a^{-1} \equiv 1 \pmod p 转化成 ax + by = c 的形式：

于是就可以使用 Exgcd 求解 a^{-1}。

::::success[Code]

int Exgcd(int a, int b, int &x, int &y) {
  if (!b) {
    x = 1, y = 0;
    return a;
  }
  ll res = Exgcd(b, a % b, y, x);
  y -= a / b * x;
  return res;
}

int main() {
  ios::sync_with_stdio(0), cin.tie(0), cout.tie(0);
  cin >> a >> b;
  Exgcd(a, b, x, y);
  cout << (x % b + b) % b << '\n';
  return 0; 
}

::::

2.2 费马小定理求逆元

• 优点：代码简单；
• 缺点：p 必须为质数；
• 时间复杂度：\mathcal{O}(\log p)。

什么是费马小定理：

如果 p 是质数，且 \gcd(a, p) = 1，那么

a^{p - 1} \equiv 1 \pmod p

现在我们开始求逆元：

由 a^{p - 1} \equiv 1 \pmod p，得：

由于 a 是常数，所以当 p 是质数时，a 的逆元为：

于是，我们可以使用快速幂计算 a^{p - 2} \bmod p。

::::success[Code]

ll Qpow(ll a, ll b, ll p) {
  ll ans = 1;
  for (; b; b >>= 1) {
    if (b & 1) {
      ans = ans * a % p;
    }
    a = a * a % p;
  }
  return ans;
}

int main() {
  ios::sync_with_stdio(0), cin.tie(0), cout.tie(0);
  cin >> n >> p;
  cout << Qpow(n, p - 2) << '\n';
  return 0; 
}

::::

2.3 线性递推求逆元

• 优点：可以一次性求 n 以内所有正整数模 p 的逆元；
• 缺点：需满足 p 为质数且 n < p；
• 时间复杂度：\mathcal{O}(n)。

令 i 的逆元为 \text{inv}_i，有如下递推式：

边界条件为：\text{inv}_1 = 1。

证明：

令 p = ki + r，其中 k = \left\lfloor \frac{p}{i} \right\rfloor，r = p \bmod i。

在模 p 意义下：

两边同乘 i^{-1}r^{-1}，得：

所以：

即：

带入 k = \left\lfloor \frac{p}{i} \right\rfloor，r = p \bmod i，得：

得证。

::::success[Code]

inv[0] = 0, inv[1] = 1;
for (int i = 2; i <= n; ++ i) {
  inv[i] = (p - p / i) * inv[p % i] % p;
} 

::::

2.4 阶乘逆元的线性求法

• 优点：可以一次性求 n 以内所有正整数的阶乘模 p 的逆元；
• 条件：p 为质数；
• 时间复杂度：\mathcal{O}(n)。

这种方法在组合数中用的比较多。

1. 先递推计算阶乘：\text{fac}_i = (\text{fac}_{i - 1} \cdot i) \bmod p（边界条件：\text{fac}_1 = 1）；
2. 用 Exgcd 或费马小定理求出 n! 的逆元 \text{invfac}_n；
3. 递推求其他阶乘逆元：\text{invfac}_i = ((i + 1) \cdot \text{invfac}_i) \bmod p。

证明：

因为 i! = \cfrac{(i + 1)!}{i + 1}，所以 (i!)^{-1} = (((i + 1)!)^{-1} \cdot (i + 1)) \bmod p。

::::success[Code]

ll Qpow(ll a, ll b) {
  ll ans = 1;
  for (; b; b >>= 1) {
    if (b & 1) {
      ans = ans * a % kMod;
    }
    a = a * a % kMod;
  }
  return ans;
}

void Init() {
  fac[0] = fac[1] = 1;
  for (int i = 2; i <= n; ++ i) {
    fac[i] = fac[i - 1] * i % kMod;
  }
  inv[n] = Qpow(fac[n], kMod - 2);
  for (int i = n - 1; i; -- i) {
    inv[i] = inv[i + 1] * (i + 1) % kMod;
  }
}

::::

Part 3：例题

其实单独考乘法逆元的题目很少，主要结合组合数学。

P1082 [NOIP 2012 提高组] 同余方程

这不是逆元板子题吗？

由于 b 未保证是质数，所以只能使用 Exgcd。

P3811 【模板】模意义下的乘法逆元

考虑使用费马小定理，但提交之后发现 TLE 了。

分析一下，费马小定理求逆元单次是 \mathcal{O}(\log p) 的，总体就是 \mathcal{O}(n \log p)，在本题数据范围下超时。

所以必须使用线性递推求逆元，时间复杂度 \mathcal{O}(n)。

P2265 路边的水沟

简单分析可得，右下角的闸门到左上角的闸门的总共有 n + m 次移动次数。

而我们必须向上走 n 次，向右走 m 次。

答案就是在 n + m 次移动中选 n 个向上的方案数，答案即为 C^{n + m}_m。

使用阶乘逆元的线性求法即可。

P5431 【模板】模意义下的乘法逆元 2

由于：

令 \text{mul}_i = \prod_{j = 1}^i a_j，则：

接下来求 a_i 的逆元：

因为：

所以：

别忘了，原式还有一个 k^i，一边算一边乘即可。

本题需要配合快读使用。

::::success[Code]

#include <bits/stdc++.h>

using namespace std;

using ll = long long;
using ull = unsigned long long;

const int kMaxN = 5e6 + 10;

ll n, p, k, a[kMaxN], mul[kMaxN], inv[kMaxN], ans;

ll R() {
  ll x = 0, f = 1;
  char ch = getchar();
  while (!isdigit(ch)) {
    if (ch == '-') {
      f = -1;
    }
    ch = getchar();
  }
  while (isdigit(ch)) {
    x = x * 10 + (ch ^ 48);
    ch = getchar();
  }
  return x * f;
}

ll Qpow(ll a, ll b) {
  ll ans = 1;
  for (; b; b >>= 1) {
    if (b & 1) {
      ans = ans * a % p;
    }
    a = a * a % p;
  }
  return ans;
}

int main() {
  ios::sync_with_stdio(0), cin.tie(0), cout.tie(0);
  n = R(), p = R(), k = R(), mul[0] = 1;
  for (int i = 1; i <= n; ++ i) {
    a[i] = R();
    mul[i] = mul[i - 1] * a[i] % p;
  }
  inv[n] = Qpow(mul[n], p - 2);
  for (int i = n - 1; i; -- i) {
    inv[i] = inv[i + 1] * a[i + 1] % p;
  }
  for (int i = n; i; -- i) {
    ans = (ans + (inv[i] * mul[i - 1]) % p) % p * k % p;
  }
  cout << ans << '\n';
  return 0;
}

::::

欧拉函数

Part 0：前置知识

• 基础数学（四则运算，符号定义等）。
• 最大公约数（gcd）；
• 互质。

Part 1：欧拉函数的定义

在信息学中，我们常常要解决互质类问题或处理大指数取模运算，而欧拉定理和欧拉函数就是一个好帮手。

定义：欧拉函数 \varphi(n) 表示小于等于 n 的正整数中与 n 互质的数的个数，即 \varphi(n) = |\{k \in \mathbb{N^+} \mid 1 \le k \le n,\ \gcd(k, n) = 1\}|。

欧拉函数有一个核心公式，此公式用处较多。

若 n = p_1^{k_1}p_2^{k_2} \cdots p_n^{k_n}，则：

展开写就是：

此公式可以使用下文的性质 1 证明。

Part 2：欧拉函数的性质

2.1：性质 1

对于任意质数 p，有 \varphi(p) = p - 1，对于 p 的幂次 p^k（k \geq 1），有 \varphi(p^k) = p^k - p^{k - 1} = p^{k - 1}(p - 1)。

证明：

在 1, 2, \cdots, p 中，不与 p 互质的数只有 p 本身，故 \varphi(p) = p - 1。

在 1, 2, \cdots, p^k 中，只有 p 的倍数不与 p 互质，即 p, 2p, 3p, \cdots, p^{k}。由于这些数有 p^{k - 1} 个，所以 \varphi(p^k) = p^k - p^{k - 1} = p^{k - 1}(p - 1)。

2.2：性质 2

若 a \mid x，则 \varphi(ax) = a \cdot \varphi(x)。

此性质证明计算量较大。

由 a \mid x，设 x = \prod_{i = 1}^n p_i^{k_i}，a = \prod_{i = 1}^n p_i^{s_i}（0 \leq s_i \leq k_i）。

那么：

使用用欧拉函数公式：

对每一项：

因此 \varphi(ax) = a \cdot \varphi(x)。

2.3：性质 3

欧拉函数是积性函数。

即，对于任意满足 \gcd(a, b) = 1 的整数 a, b，都有 \varphi(ab) = \varphi(a)\varphi(b)。

特别地，当 n \bmod 2 = 1 时，\varphi(2n) = \varphi(n)。

此性质证明需要用到作者的知识盲区，作者太弱所以就不写上去了。

Part 3：欧拉函数的计算

3.1：暴力计算欧拉函数

欧拉函数公式，枚举每一个质因子然后计算答案即可。

时间复杂度：单个 \mathcal{O}(\log n)。

::::success[Code]

ll Phi(ll n) {
  ll ans = n;
  for (ll i = 2; i * i <= n; ++ i) {
    if (!(n % i)) {
      ans -= ans / i;
      for (; !(n % i); n /= i);
    }
  }
  if (n > 1) {
    ans -= ans / n;
  }
  return ans;
}

::::

3.2 线性筛求解欧拉函数

我们借助线性筛来求解欧拉函数，这个求法主要依赖于欧拉函数的性质。先放代码，再讲原理。

::::success[Code]

const int kMaxN = 1e6 + 10;

ll p[kMaxN], phi[kMaxN], cnt;
bool vis[kMaxN];

void Init() {
  vis[1] = phi[1] = 1;
  for (ll i = 2; i < kMaxN; ++ i) {
    if (!vis[i]) {
      phi[i] = i - 1;
      p[++ cnt] = i;
    }
    for (ll j = 1; j <= cnt && i * p[j] < kMaxN; ++ j) {
      vis[i * p[j]] = 1;
      phi[i * p[j]] = phi[i] * (p[j] - 1);
      if (!(i % p[j])) {
        phi[i * p[j]] = phi[i] * p[j];
        break;
      }
    }
  }
}

::::

令 p_i 表示第 i 个质数，\text{vis}_i = 1 表示 i 是合数。

考虑枚举 2 至 N 的所有整数。

如果 \text{vis}_i = 1，就证明 i 未被筛出，即 i 是质数，那么 \varphi(i) = i - 1。

然后，我们就开始筛除 i 的倍数。如果 i \mid p_j，那么根据性质 2，得到 \varphi(i \cdot p_j) = \varphi(i) \cdot p_j。否则，我们根据性质 3（积性函数的性质），可以得到 \varphi(i \cdot p_j) = \varphi(i) \cdot (p_j - 1)。

时间复杂度 \mathcal{O}(N)。

Part 4：例题

UVA11327 Enumerating Rational Numbers

解读一下代码，就可以发现对于每一个 d，满足条件的分数个数就是 \varphi(d)。

于是尝试枚举 d，每一次将个数与 \varphi(d) 相加。如果发现超过 k 了，那就证明答案的 d 一定是当前的 d，此时枚举 n 即可。

由于 d \le 2 \times 10^5（由样例发现），线性筛预处理即可。

CF1295D Same GCDs

首先对式子做一次辗转相除，得：

显然题目可以转化为求有多少 x 满足 \gcd(x, m) = \gcd(a, m)。

令 \gcd(x, m) = \gcd(a, m) = d，则 \gcd\left(\cfrac{x}{d}, \cfrac{m}{d}\right) = 1。

现在再次转换题目：求与 \cfrac{m}{d} 互质的数的个数。显然就是欧拉函数，答案为 \varphi \left(\cfrac{m}{d}\right)。

P1891 疯狂 LCM

遇到这种题，无脑暴力推式子：

分析上式，目前的目标就是求 \sum_{i = 1}^{d} i \cdot [\gcd(i, d) = 1]。

假设你欧拉函数掌握的不是特别好，我们来找一下规律。假设 n = 8，那么 \gcd(i, n) = 8 的有 1, 3, 5, 7。这时细心的人可以发现 1 + 7 = 3 + 5 = n = 8。

没错！对于 \gcd(i, d) = 1，也存在 \gcd(d - i,\ d) = 1（证明很简单）！

那么：

答案即为：

那么，我们先 \mathcal{O}(N) 预处理欧拉函数，然后 \mathcal{O}(N \log N) 预处理每一个 d 的 f(d) = \sum_{d \mid n} d \cdot \frac{\varphi(d)}{2}（可以使用 DP 处理）。最后，询问时直接输出 n \cdot f(n) 即可。

::::success[Code]

#include <bits/stdc++.h>

using namespace std;

using ll = long long;
using ull = unsigned long long;

const int kMaxN = 1e6 + 10;

ll n, phi[kMaxN], p[kMaxN], f[kMaxN], cnt;
bool vis[kMaxN];

void Init() {
  vis[1] = phi[1] = 1;
  for (ll i = 2; i <= 1e6; ++ i) {
    if (!vis[i]) {
      phi[i] = i - 1;
      p[++ cnt] = i;
    }
    for (ll j = 1; j <= cnt && i * p[j] <= 1e6; ++ j) {
      vis[i * p[j]] = 1;
      phi[i * p[j]] = phi[i] * (p[j] - 1);
      if (!(i % p[j])) {
        phi[i * p[j]] = phi[i] * p[j];
        break;
      }
    }
  }
  for (ll i = 1; i <= 1e6; ++ i) {
    for (ll j = i; j <= 1e6; j += i) {
      f[j] += (phi[i] * i + 1) / 2;
    }
  }
}

int main() {
  ios::sync_with_stdio(0), cin.tie(0), cout.tie(0);
  Init();
  int t;
  for (cin >> t; t; -- t) {
    cin >> n;
    cout << n * f[n] << '\n';
  }
  return 0;
}

::::

中国剩余定理（CRT）

Part 0：前置知识

• 基础数学（四则运算，符号定义等）。
• 互质；
• 同余；
• 逆元；
• 裴蜀定理。

Part 1：中国剩余定理是什么？

中国剩余定理（Chinese Remainder Theorem，CRT），又称孙子定理，是数论中最具代表性的定理之一，最早记载于我国南北朝时期的数学著作《孙子算经》，其中经典的“物不知数”问题便是其雏形：“今有物不知其数，三三数之剩二，五五数之剩三，七七数之剩二，问物几何？”

中国剩余定理一般用来求解以下一次线性同余方程组：

其中对于 i \neq j，有 \gcd(m_i, m_j) = 1。

Part 2：同余方程组的求法

我们考虑使用余数的可加性，那么题目转化为求出 n 个方程组的解（i \in [1, n]）：

那么 x = \sum_{i = 1}^n x_i。

进一步转化，其实就是求这 n 个方程组的解：

由余数的可乘性，得 x_i = y_ib_i。由于对于 j \neq i，有 y_i \equiv 0 \pmod{m_j}，那么 y_i 一定可以表示为 k_i \cdot \prod_{j = 1}^n m_i[j \neq i]（k \in \mathbb{N^+}）。

又 y_i \equiv 0 \pmod{m_1}，那么 k_i \cdot \prod_{j = 1}^n m_i[j \neq i] \equiv 1 \pmod {m_i}。观察式子，发现与逆元的定义相同。那么求 y_i 就转为了求 \prod_{j = 1}^n m_j[j \neq i] \bmod m_i 的逆元。这时，就可以拿 Exgcd 求解了。

求出 n 个方程组的解后，就得到了原方程的整数解。但因为题目要求最小非负整数解，所以还要再进行对乘积的取模。

::::success[Code]

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 CRT() {
  ll mul = 1, ans = 0;
  for (int i = 1; i <= n; ++ i) {
    mul *= b[i];
  }
  for (ll i = 1, k, x, y; i <= n; ++ i) {
    k = mul / b[i];
    Exgcd(k, b[i], x, y);
    ans = (ans + k * b[i] * x % mul) % mul;
  }
  return (ans % mul + mul) % mul;
}

::::

Part 3：扩展中国剩余定理（ExCRT）

由于中国剩余定理的局限性太强（两两互质），所以我们需要扩展中国剩余定理。

在信息学中，CRT 与 ExCRT 复杂度一致。所以，只要学习扩展中国剩余定理就可以覆盖 CRT 的问题。

考虑如下方程组：

其中 m_i 是任意正整数。

首先，我们将式子等价转换。

对于 x \equiv b_1 \pmod{m_1} 和 x \equiv b_2 \pmod{m_2}，则有：

于是：

移项，得：

观察此式子，发现形态与线性不定方程相同，考虑使用扩展欧几里得（Exgcd）求解。

令 d = \gcd(m_1, m_2)。

首先，根据裴蜀定理：若 d \nmid (b_2 - b_1)，方程无解。否则，我们可以使用 Exgcd 求解。

设通过 Exgcd 得到：

两边同乘 \frac{b_2 - b_1}{d}：

得到一组特解：

那么，原方程组的通解为：

然后我们考虑使用数学归纳法。

设前 k - 1 个方程的一个特解为 x'，则通解为 x'+ t \cdot \text{lcm}(m_1, m_2, \cdots, m_{k - 1})（t \in \mathbb{Z^+}）。对于第 k 个方程 x \equiv b_k \pmod{m_k}，将通解代入第 k 个方程的 x：

移项，得：

将同余方程转成一次不定方程，得：

若该方程无解，则说明前 k 个方程无解，退出。否则，更新特解 x_1 和 \text{lcm}。

::::success[Code]

ll Exgcd(ll a, ll b, ll &x, ll &y) {
  if (!b) {
    x = 1, y = 0;
    return a;
  }
  ll res = Exgcd(b, a % b, y, x);
  y -= a / b * x;
  return res;
}

ll Qmul(ll a, ll b, ll p) {
  ll ans = 0;
  for (; b; b >>= 1) {
    if (b & 1) {
      ans = (ans + a) % p;
    }
    a = (a + a) % p;
  }
  return ans;
}

ll ExCRT() {
  ans = b[1], lcm = m[1];
  for (ll i = 2, d, x, y, k; i <= n; ++ i) {
    b[i] = ((b[i] - ans) % m[i] + m[i]) % m[i];
    d = Exgcd(lcm, m[i], x, y);
    if (b[i] % d) {
      ans = -1;
      break;
    }
    k = Qmul(x, b[i] / d, m[i]);
    ans += k * lcm, lcm = lcm / d * m[i];
    ans = (ans % lcm + lcm) % lcm;
  }
  return ans;
}

::::

Part 4：例题

P1495 【模板】中国剩余定理（CRT）/ P4777 【模板】扩展中国剩余定理（EXCRT）

板子题，套模板即可。

注意 P1495 需要 __int128。

P3868 [TJOI2009] 猜数字

题目中的式子很难看，那我们就把它写成同余方程组：n - a_i \equiv 0 \pmod{b_i}。

因为 a \equiv b \pmod m 时，有 a + c \equiv b + c \pmod m。所以题目中的式子就转化成 n \equiv a_i \pmod{b_i}。

然后就是模板了。

CF687B Remainders Game

此题思维性较强。

我们由题目可以知道，Arya 可以从 Pari 手中得知 x \bmod c_i 的值。也就是，令第 i 次询问的结果是 b_i，那么 Alice 可以得到如下式子：

这时你会发现，这个方程组与 ExCRT 一摸一样！那么，我们结合 ExCRT，可以发现，通过任意的：

可以得到 x \bmod \text{lcm}(c_i, c_{i + 1}) 的值！也就是，我们只需要知道 c_i, c_{i + 1} 和 x \bmod c_i，x \bmod c_{i + 1}，就可以得出 x \bmod \text{lcm}(c_i, c_{i + 1}) 的值。

因此，我们只需要求出 \text{lcm}(c_1, c_2, \cdots, c_n) 的值，然后判断是否有 k \mid \text{lcm}(c_1, c_2, \cdots, c_n) 即可。

为什么呢？因为当 \text{lcm}(c_1, c_2, \cdots, c_n) 满足上述条件时，我们可以得到 x \bmod tk（t 为任意正整数）的值。那么，根据同余的性质，便可以得到 x \bmod k。

欧拉定理

Part 0：前置知识

• 基础数学（四则运算，符号定义等）。
• 互质。

Part 1：欧拉定理的内容

1.1：欧拉定理

若正整数 a, n 满足 \gcd(a, n) = 1，那么有：

当 n 为质数 p 时，\varphi(p) = p - 1，欧拉定理退化为费马小定理：

1.2：证明

接下来，我们证明欧拉定理。

但是，在证明之前，我们需要先了解一个关键的观察：

取 n = 8，a = 3，将与 n 互质的数乘以 a 在对 n 取模：

我们发现，得到的余数是 3, 1, 7, 5，恰好是与 n 互质的数 1, 3, 5, 7 的一个重新排列。

这不是巧合，这是一个普遍规律。

那么，接下来我们开始证明。

设所有小于等于 n 且与 n 互质的正整数为：

由欧拉函数，得 k = \varphi(n)。

考虑将这些数分别乘上 a，然后对 n 取模：

接下来，我们需要证明三个性质。

性质 1：\gcd(ar_i, n) = 1。

因为 \gcd(a, n) = \gcd(r_i, n) = 1，那么根据互质的性质，得 \gcd(ar_i, n) = 1。

性质 2：对于 i \neq j，有 ar_i \not \equiv ar_j \pmod n。

假设 ar_i \equiv ar_j \pmod n，即 n \mid a(r_i - r_j)。

因为 \gcd(a, n) = 1，所以 n \mid (r_i - r_j)，即 r_i \equiv r_j \pmod n。

又 r_i, r_j \in [1, n] 且 \gcd(r_i, n) = \gcd(r_j, n) = 1，所以 r_i = r_j。

因为 r_i, r_j 互不相同，所以 i = j。

与 i \neq j 矛盾，原命题得证。

性质 3：(ar_1 \bmod n) \in [1, n)。

将命题转化成证明 ar_1 \bmod n \neq 0，即 ar_1 \not \equiv 0 \pmod n。

假设 ar_1 \equiv 0 \pmod n，那么 n \mid ar_i，但 \gcd(a, n) = \gcd(r_i, n) = 1。

矛盾，原命题得证。

三个性质证明完后，我们会发现一个关键结论。

现在我们有两个集合：

• 原始集合 S = \{r_1, r_2, \cdots, r_{\varphi(n)}\}；
• 乘积集合 S' = \{ar_1 \bmod n,\ ar_2 \bmod n,\ \cdots,\ ar_{\varphi(n)} \bmod n\}。

因为乘积集合也具有以下三个性质：

• 对于 i \neq j，有 ar_i \not \equiv ar_j \pmod n；

这表明一个关键结论：S = S'。

既然 S = S'，那么元素的乘积应该相等，那么推导出：

将 \prod_{i = 1}^{\varphi(n)} ar_i \bmod n 展开：

所以有：

由于原始集合的乘积与 n 互质，故原始集合的乘积在模 n 下有乘法逆元，在式子两边同乘它的逆元：

证毕。

Part 2：欧拉定理的应用

2.1：求乘法逆元

如果 \gcd(a, n) = 1，则 a 在模 n 意义下的乘法逆元为：

证明：

::::warning[注意] 此处实际上比费马小定理求逆元的应用更广。 ::::

2.2：指数降幂

对于 a^b \bmod n，如果 \gcd(a, n) = 1，则：

证明：

令 b = q \cdot \varphi(n) + r（其中 r \in [0, \varphi(n))），则：

Part 3：扩展欧拉定理

对于任意正整数 n 和任意非负整数 b，有：

注意：不需要满足 \gcd(a, n) = 1。

扩展欧拉定理给出了降幂更加广泛的形式。

::::success[Code]

#include <bits/stdc++.h>

using namespace std;

using ll = long long;
using ull = unsigned long long;

const int kMaxN = 1e6 + 10;

ll a, b, m, kMod;
bool ok;

ll R() {
  ll k = 0, f = 1;
  char c = getchar();
  while (c < '0' || c > '9') {
    (c == '-') && (f = -1);
    c = getchar();
  }
  while (c >= '0' && c <= '9') {
    k = k * 10 + c - '0';
    if (k >= kMod) {
      ok = 1;
    }
    k %= kMod;
    c = getchar();
  }
  return k * f;
}

ll Phi(ll n) {
  ll ans = n;
  for (ll i = 2; i * i <= n; ++ i) {
    if (!(n % i)) {
      ans -= ans / i;
      for (; !(n % i); n /= i);
    }
  }
  if (n > 1) {
    ans -= ans / n;
  }
  return ans;
}

ll Qpow(ll b, ll p) {
  ll ans = 1;
  for (; p; p >>= 1) {
    if (p & 1) {
      ans = ans * b % m;
    }
    b = b * b % m;
  }
  return ans;
}

int main() {
  cin >> a >> m;
  kMod = Phi(m), b = R();
  if (ok) {
    b += kMod;
  }
  cout << Qpow(a, b) << '\n';
  return 0;
}

::::

Part 4：例题

P5091 【模板】扩展欧拉定理

板子题，套模板即可。

P4139 上帝与集合的正确用法

上面的文字很多，不用在意。我们直接将这个数列转化成求：

也就是 2 的无限幂塔。

考虑使用扩展欧拉定理。可以证明，指数满足 b \ge varphi(p)。

于是，考虑使用递归的形式，直到边界条件 p = 1，此时式子的值是 0。

P2350 [HAOI2012] 外星人

由于一般的题目不单独考欧拉定理，所以此题严格上归于“欧拉函数”。

题目的意思就是：给定 n，每次将 n 变为 \varphi(n)，求几次操作后 n = 1。

首先，人人皆知，只有 n = 1 或 2 时才有 \varphi(n) = 1。

接下来，根据题目中给的提示，我们有：

观察到，式子中有一个 (p_i - 1)。那么，每一次会把上一次操作的答案中的最多一个 2 变为 1。

所以，我们只需要求出 n 的每一个质因子会产生多少个 2 即可。令 n = i 时答案为 f_i，我们分情况讨论：

• 当 i 是质数时，显然有递推式 f_i = f_{i - 1}；
• 当 i 不是质数时，令 i = p \cdot q。那么根据乘积的性质，可以得到 f_i = f_p + f_q。

你以为这就结束了？错！如果当原数 n 中没有因子 2 时，第一次操作不会把任何 2 变为 1。

所以，当原数 n 中没有因子 2 时，答案需减少 1。

::::success[Code]

#include <bits/stdc++.h>

using namespace std;

using ll = long long;
using ull = unsigned long long;

const int kMaxN = 1e5 + 10;

ll m, f[kMaxN], p[kMaxN], q[kMaxN], cnt, ans;
bool isp[kMaxN];

void Init() {
  f[1] = 1;
  for (ll i = 2; i < kMaxN; ++ i) {
    if (!isp[i]) {
      p[++ cnt] = i, f[i] = f[i - 1];
    }
    for (ll j = 1; j <= cnt && p[j] * i < kMaxN; ++ j) {
      isp[p[j] * i] = 1, f[p[j] * i] = f[p[j]] + f[i];
      if (!(i % p[j])) {
        break;
      }
    }
  }
}

int main() {
  ios::sync_with_stdio(0), cin.tie(0), cout.tie(0);
  int t;
  Init();
  for (cin >> t; t; -- t) {
    cin >> m, ans = 1;
    for (int i = 1; i <= m; ++ i) {
      cin >> p[i] >> q[i];
      ans += f[p[i]] * q[i] + (p[i] % 2? 0 : -1);
    }
    cout << ans << '\n';
  }
  return 0;
}

::::

小步大步算法（BSGS）

Part 0：前置知识

• 基础数学（四则运算，符号定义等）。
• 模运算；
• 欧拉定理；
• 离散对数（可以不掌握）。

Part 1：小步大步算法

1.1：算法解决的问题

小步大步算法（Baby Step Giant Step）一般解决如下问题：

给定一个质数 p 和整数 a, b（其中 a \in (0, p)，b \in [0, p)），求最小的非负整数 x，满足：

或者报告无解。

1.2：算法思想

此算法的思想是根号类思想，结合 Meet in middle 的思想。

令 t = \lceil\sqrt{p}\rceil，用 t 表示 x，则 x = At - B，那么：

两边同乘 a^B：

考虑枚举 B，范围是 [0, t - 1]，将所有的 b \cdot a^{B} 用 Hash 表存储。

再考虑枚举 A，范围是 [1, t]，计算 A 个 a^t 的积。如果在 Hash 表中出现过，说明找到了对应的 A, B。

此时，x = At - B 即为答案。

总时间复杂度 \mathcal{O}(t)，即 \mathcal{O}(\sqrt{p})。

::::success[Code]

#include <bits/stdc++.h>

using namespace std;

using ll = long long;
using ull = unsigned long long;

const int kMaxN = 1e5 + 10;

ll p, b, n;

ll Qpow(ll a, ll b, ll p) {
  ll ans = 1;
  for (; b; b >>= 1) {
    if (b & 1) {
      ans = ans * a % p;
    }
    a = a * a % p;
  }
  return ans;
}

ll BSGS(ll p, ll b, ll n) {
  unordered_map<ll, ll> mp;
  ll t = sqrt(p) + 1;
  for (ll i = 1, mul = (n * b) % p; i <= t; ++ i, mul = (mul * b) % p) {
    mp[mul % p] = i;
  } 
  ll pw = Qpow(b, t, p);
  for (ll i = 1, mul = pw; i <= t; ++ i, mul = (mul * pw) % p) {
    if (mp[mul]) {
      return t * i - mp[mul];
    }
  }
  return -1;
} 

int main() {
  ios::sync_with_stdio(0), cin.tie(0), cout.tie(0);
  cin >> p >> b >> n;
  ll ans = BSGS(p, b, n);
  if (ans == -1) {
    cout << "no solution\n";
  } else {
    cout << ans << '\n';
  }
  return 0;
}

::::

Part 2：扩展小步大步算法

在扩展小步大步算法（exBSGS）中，不要求 \gcd(a, p) = 1。

那么我们开始推导。

首先，我们尝试将同余方程写成不定方程：

根据裴蜀定理，当 \gcd(a, p) \nmid b 时方程无整数解。否则，我们令 d = \gcd(a, b)，将不定方程两边同除 d：

再将不定方程写成同余方程的形式：

两边同乘 \frac{a}{d} 在模 \frac{p}{d} 下的逆元：

我们发现，这跟 BSGS 的形式完全相同。

所以，我们检查 \gcd(a, \frac{p}{d}) 是否为 1。如果满足，则使用 BSGS 求出 x - 1。否则，重复执行上述操作，直到满足条件或得到无解。

::::success[Code]

#include <bits/stdc++.h>

using namespace std;

using ll = long long;
using ull = unsigned long long;

const int kMaxN = 2e5 + 10;

ll BSGS(ll a, ll p, ll b, ll pw) {
  unordered_map<ll, ll> mp;
  ll t = sqrt(p) + 1, mul = 1;
  for (ll i = 1; i <= t; ++ i, mul = mul * a % p) {
    mp[mul * b % p] = i;
  }
  for (ll i = 1, k = mul, mul = pw; i <= t + 1; ++ i, mul = mul * k % p) {
    if (mp[mul] && t * (i - 1) - mp[mul] + 1 >= 0) {
      return t * (i - 1) - mp[mul] + 1;
    }
  }
  return -1;
}

ll exBSGS(ll a, ll p, ll b) {
  a %= p, b %= p;
  if (b == 1 || p == 1) {
    return 0;
  }
  ll res = 0, mul = 1;
  for (ll d; (d = __gcd(a, p)) != 1; ) {
    if (b % d) {
      return -1;
    }
    ++ res, b /= d, p /= d;
    mul = (mul * a / d) % p;
    if (mul == b) {
      return res;
    } 
  }
  ll ans = BSGS(a, p, b, mul);
  return (ans == -1? ans : ans + res);
}

int main() {
  ios::sync_with_stdio(0), cin.tie(0), cout.tie(0);
  for (ll a, p, b, ans; ; ) {
    cin >> a >> p >> b;
    if (!a) {
      return 0;
    }
    ans = exBSGS(a, p, b);
    if (ans == -1) {
      cout << "No Solution\n";
    } else {
      cout << ans << '\n';
    }
  }
  return 0; 
}

::::

Part 3：例题

P3846 【模板】BSGS / P4195 【模板】扩展 BSGS

模板题，套板子即可。

P4454 [CQOI2018] 破解D-H协议

首先，将题目中 2、3 步的式子写成同余方程形式：

这不就是 BSGS 的形式了？于是使用 BSGS 求出 a, b 之后计算 g^{ab} \bmod P 即可。

P4884 多少个 1？

感觉这一题难度虚高。

首先，对 N 个 1 的数推导通项公式：

那么方程转化为：

两边同乘 9，得：

移项，得：

于是就转化为标准 BSGS 形式了。