题解 - ZR 2590 朝花夕拾 / 模质数意义下质数次高次剩余

题意简述

给出自然数 $a$ 和质数 $p,q$，求解 $x^q\equiv a\pmod p$，输出任一符合条件的 $x\in [0,p)$ 即可，无解输出 -1。

$T\le 100$ 组数据，$p,q\le 10^{18}$，$0\le a<p$。

解法

接下来的题解会尽量依照正向的思路去写，笔者就是被反向思路写的官方题解打蒙的……

本题构造思路太神仙了，笔者认为实用性不强，在考场上可能很难自主复现，估计也考不着。Anyway, % 127。

简单解

出题人说这是高次剩余板子，可是看数据范围，朴素离散对数显然不能过，时间复杂度达到 $O(T\sqrt p)$。如果我们需要上 BSGS 算法，模数不可以是 $p$ 级别的。

考虑费马小定理，有 $a\equiv a^{m(p-1)+1}\pmod p$。

那么就可能有简单解：如果存在 $m$ 满足 $q|m(p-1)+1$，则方程有解

$$x\equiv a^{\frac{m(p-1)+1}{q}} \pmod p$$

就是直接在指数上除以 $q$。

考虑出现这种情况的条件：发现 $q\perp p-1$，上述解即化为

$$x\equiv a^{q^{-1}\pmod{p-1}} \pmod p$$

判无解

然后判断一下是否无解。

为了便于判断，令 $q^Kt=p-1$（显然指数应该放在质数 $q$ 上）。判断完上面 $K=0$ 的情况，接下来我们有 $K\ge 1$，$q|(p-1)$。

若有解 $x$，则

$$x^{p-1}\equiv 1\pmod p$$

$$a^{\frac{p-1}{q}}\equiv 1\pmod p$$

$$a^{q^{K-1}t}\equiv 1\pmod p$$

若不满足上述条件，则无解。

由判定的同余式得到的启示

发现上面的写法和费马小定理的形式有些像，考虑用类似于简单解的形式，找出通用解。

假设我们已经有 $k$ 满足 $a^{q^kt}\equiv 1\pmod p$（注意区分 $K$ 和 $k$），且能找到 $m$ 满足 $q|mq^kt+1$，则由 $a\equiv a^{mq^kt+1}\equiv (a^{\frac{mq^kt+1}{q}})^q\pmod p$，有解

$$x\equiv a^{\frac{mq^kt+1}{q}}\pmod p$$

发现能找到 $m$ 的条件是 $q^kt\perp q$。所以 $k=0$，解化为

$$x\equiv a^{q^{-1}\pmod{t}}\pmod p$$

由 $a^{q^{K-1}t}\equiv 1\pmod p$，我们已经有一个满足假设的 $k=K-1$。我们只要让它不断地减小，直到它等于 $0$。在过程中，$a^{q^{K-1}t}$ 最后会变化为 $a_0$，最终有 $a_0^t\equiv 1\pmod p$，容易求解 $x_0^t\equiv a_0\pmod p$。求出最后的解 $x_0$，我们再以一种方式回带即可。

当然，在 $K=1$ 时我们就不需要再去减小了，因此接下来 $K\ge 2$，即 $q\le\sqrt p$。

（讲完算法之后，会发现这个位置对时间复杂度分析是重要的，因为这时 $q$ 已经较小，我们可以直接对满足 $w^q\equiv 1\pmod p$ 的 $w$ 做 BSGS。）

让 $k$ 减一

我们已经有 $a^{q^kt}\equiv 1\pmod p$。令 $v\equiv a^{q^{k-1}t}\pmod p$，则 $v^q\equiv a^{q^kt}\equiv 1\pmod p$。

你可以判一下，如果此时的 $v=1$，那么 $k$ 减一时不用做其他变化，直接 continue 即可。

否则我们希望找到一个 $a'^{q^{k-1}t}\equiv a^{q^kt}\equiv 1\pmod p$，用 $a'$ 替换 $a$。

考虑用 $a'\equiv z^qa\pmod p$ 来描述这种变换（这样子设，回带 $x$ 是简明的）。

$$a'^{q^{k-1}t}\equiv (z^qa)^{q^{k-1}t}\equiv z^{q^kt}a^{q^{k-1}t}\equiv z^{q^kt}v\equiv 1\pmod p$$

$$z^{q^kt}\equiv v^{-1}\pmod p$$

由于 $q$ 小，我们想能不能直接找到一个 $Z$，对它求 $b$ 满足 $Z^b\equiv v^{-1}\pmod p$（假设一），并且 $Z^q\equiv 1\pmod p$（假设二）……？

在构造 $Z$ 的路上狂奔

回到最初的费马小定理，对于任一个 $w$，我们都有 $w^{q^Kt}\equiv w^{p-1}\equiv 1\pmod p$。要满足假设二，我们取一个

$$Z\equiv w^{q^{K-1}t}\pmod p$$

即可。

假设一我们尝试用 BSGS 算法求出 $b$ 来完成。由 $Z^q\equiv 1\pmod p$ 联系到 $v^q\equiv 1\pmod p$ 这一重要的性质，可以证明 $Z\neq 1$ 且 $v\neq 1$ 时都可以求出 $b$。也就是说，我们在选取 $w$ 时，要令 $w^{q^{K-1}t}\not\equiv 1\pmod p$。

然后我们来用 $Z$ 生成 $z$。

$$z^{q^kt}\equiv Z^b\equiv w^{q^{K-1}tb}\pmod p$$

$$z\equiv w^{q^{K-1-k}b}\pmod p$$

然后再用 $z$ 计算 $a'$。如此循环，让 $k$ 不断减一，直到 $0$。

回带

我们最终求出了 $x_0^q\equiv a_0\pmod p$ 的一个解，我们要回去求 $x^q\equiv a\pmod p$。

根据 $a_0\equiv z_0^qa_1\pmod p$，我们有

$$x_0^q\equiv z_0^qa_1\pmod p$$

$$(\frac{x_0}{z_0})^q\equiv a_1\pmod p$$

所以可以有解

$$x_1\equiv\frac{x_0}{z_0}\pmod p$$

依此类推，可以有

$$x\equiv\frac{x_0}{\prod z}\pmod p$$

时间复杂度分析

对于特判等细节处，显然是 $O(\log p)$ 的。

然后我们做的第一件非 $O(\log p)$ 的事是找到一个合法的 $w$。回顾对 $w$ 的限制：$w^{q^{K-1}t}\equiv w^{\frac{p-1}{q}}\not\equiv 1\pmod p$。显然质数 $p$ 的原根就满足这个性质，原根的指数小于 $p-1$ 的话 $\bmod p$ 后必不等于 $1$。所以我们找到的 $w$ 完全可以不超过 $p$ 的最小原根，大小是 $O(p^{\frac{1}{4}})$ 的。

接下来，我们不断让 $k$ 减一，显然 $k$ 的级别是 $O(\log_q p)$ 的。然后，我们要跑若干次 BSGS。注意到次数小于 $q$，并且要做 BSGS 的条件是 $K\ge 2$，$q\le\sqrt p$。所以 BSGS 单次复杂度 $O(\sqrt q)$ 也是 $O(p^{\frac{1}{4}})$。

综上，时间复杂度是 $O(Tp^{\frac{1}{4}}\log p)$ 的。看起来好像不能光速幂，而且这个复杂度已经可以过了。

代码

#include <bits/stdc++.h>
using namespace std;

typedef long long ll;
typedef __int128_t i128;
int T;
ll A, Q, P;
ll gcd(ll x, ll y) {
  while (y) {
    ll t(x % y);
    x = y;
    y = t;
  }
  return x;
}
ll pw(ll b, ll p, ll m) {
  ll a(1);
  while (p) {
    if (p & 1) a = static_cast<i128>(a) * b % m;
    b = static_cast<i128>(b) * b % m;
    p >>= 1;
  }
  return a;
}
void exgcd(ll a, ll b, ll& x, ll& y) {
  if (!b) x = 1, y = 0;
  else exgcd(b, a % b, y, x), y -= a / b * x;
}
ll inv(ll a, ll p) {
  ll x(0), y(0);
  exgcd(a, p, x, y);
  return (x % p + p) % p;
}
class BSGSA {
 private:
  unordered_map<ll, int> mp;
  ll bse, dvs, rng;
 public:
  BSGSA(ll b, ll p, ll m) : bse(b), dvs(m), rng(p) {
    ll a(pw(b, rng, m));
    for (ll i(1), j(a); i <= rng; ++i, j = static_cast<i128>(j) * a % m)
      mp[j] = i;
  }
  int query(ll x) {
    for (ll i(0), j(x);; ++i, j = static_cast<i128>(j) * bse % dvs)
      if (mp[j])
        return mp[j] * rng - i;
    return -1;
  }
};
ll solve(ll a, ll q, ll p) {
  ll K(0), t(p - 1);
  while (t % q == 0) ++K, t /= q;
  if (K == 0) return pw(a, inv(q, p - 1), p);
  if (pw(a, (p - 1) / q, p) != 1) return -1;
  ll iv(1);
  if (K > 1) {
    ll Z, w(1);
    while ((Z = pw(w, (p - 1) / q, p)) == 1) ++w;
    BSGSA bsgs(Z, sqrt(q) + 2, p);
    for (int k(K - 1); k; --k) {
      ll v(pw(a, t * pw(q, k - 1, p), p));
      if (v == 1) continue;
      ll b(bsgs.query(inv(v, p)));
      ll z(pw(w, b * pw(q, K - k - 1, p), p));
      a = static_cast<i128>(a) * pw(z, q, p) % p;
      iv = static_cast<i128>(iv) * z % p;
    }
  }
  return static_cast<i128>(pw(a, inv(q, t), p)) * inv(iv, p) % p;
}
int main() {
  ios::sync_with_stdio(false);
  cin.tie(nullptr);
  cin >> T;
  while (T--) {
    cin >> A >> Q >> P;
    cout << solve(A, Q, P) << endl;
  }
  return 0;
}

---

要特判 $a=0$，输出 $0$。

Updated on 2023.8.22.