初等数论（二）

数论分块

数论分块可以快速计算一些含有除法向下取整的和式（即形如 $\sum_{i=1}^nf(i)g(\left\lfloor\dfrac ni\right\rfloor)$ 的和式）。当可以在 $O(1)$ 内计算 $f(r)-f(l)$ 或已经预处理出 $f$ 的前缀和时，数论分块就可以在 $O(\sqrt n)$ 的时间内计算上述和式的值。

long long H(int n) {
  long long res = 0;  
  int l = 1, r;      
  while (l <= n) {
    r = n / (n / l); //与n/l相等的最大的r
    res += (r - l + 1) * 1LL * (n / l);  
    l = r + 1;  
  }
  return res;
}

莫比乌斯反演

设 $f(n),g(n)$ 为两个数论函数。

形式一：如果有 $f(n)=\sum_{d\mid n}g(d)$，那么有 $g(n)=\sum_{d\mid n}\mu(d)f(\frac{n}{d})$。

形式二：如果有 $f(n)=\sum_{n|d}g(d)$，那么有 $g(n)=\sum_{n|d}\mu(\frac{d}{n})f(d)$。

其中，$\mu$ 为莫比乌斯函数，定义为

$$\mu(n)= \begin{cases} 1&n=1\\ 0&n\text{ 含有平方因子}\\ (-1)^k&k\text{ 为 }n\text{ 的本质不同质因子个数}\\ \end{cases}$$

莫比乌斯函数不仅是积性函数，还有如下性质：

$$\sum_{d\mid n}\mu(d)= \begin{cases} 1&n=1\\ 0&n\neq 1\\ \end{cases}$$

即 $\displaystyle\sum_{d\mid n}\mu(d)=[n=1]$，可推出$\displaystyle [\gcd(i,j)=1]=\sum_{d\mid\gcd(i,j)}\mu(d)$

更多：莫比乌斯反演-让我们从基础开始

线性筛莫比乌斯函数

int mu[maxn], lp[maxn], p[maxn], tot;
void getmu() {
    mu[1] = 1;
    for (int i = 2;i < maxn;i++) {
        if (!lp[i]) lp[i] = p[++tot] = i, mu[i] = -1;
        for (int j = 1;j <= tot && p[j] <= lp[i] && p[j] * i < maxn;++j) {
            lp[p[j] * i] = p[j];
            if (i % p[j] == 0) {
                mu[i * p[j]] = 0;
                break;
            }
            mu[i * p[j]] = -mu[i];
        }
    }
}

中国剩余定理

中国剩余定理可求解如下形式的一元线性同余方程组（其中 $n_1, n_2, \cdots, n_k$ 两两互质）：

$$\begin{cases} x &\equiv a_1 \pmod {n_1} \\ x &\equiv a_2 \pmod {n_2} \\ &\vdots \\ x &\equiv a_k \pmod {n_k} \\ \end{cases}$$

理解：

可以设：$x \equiv e_1a_1+e_2a_2+...e_ka_k \pmod {M}$

其中, $M=n_1*n_2*...n_k$，则$e_k\equiv 0\pmod{M/n_k}$，为使上式满足方程里的所有条件，还可以推出$e_k\equiv 1\pmod{n_k}$，用扩展欧几里得算法可求得$e_k$。

void exgcd(int a, int b, int& d, int& x, int& y) {
	if (!b) { d = a; x = 1; y = 0; }
	else { exgcd(b, a % b, d, y, x); y -= x * (a / b); }
}
int mul_mod(int a, int b, int p) {
	if (b == 0)
		return 0;
	if (b == 1)
		return a % p;
	int temp = mul_mod(a, b >> 1, p);
	if (b & 1)
		return ((temp + temp) % p + a) % p;
	return (temp + temp) % p;
}
int china(int a[], int b[], int len)
{
	int ans = 0, lcm = 1, d, x, y;
	for (int i = 1;i <= len;++i) {
		lcm *= b[i];
		a[i] = (a[i] % b[i] + b[i]) % b[i];
	}
	for (int i = 1;i <= len;++i)
	{
		int tp = lcm / b[i];
		exgcd(tp, b[i], d, x, y);
		x = (x % b[i] + b[i]) % b[i];
		ans = (ans + mul_mod(mul_mod(tp, x, lcm), a[i], lcm)) % lcm;
	}
	return (ans + lcm) % lcm;
}

扩展中国剩余定理

当模数不两两互质时，可用扩展剩余定理求解。

设两个方程分别是 $x\equiv a_1 \pmod {m_1}$、$x\equiv a_2 \pmod {m_2}$；

将它们转化为不定方程：$x=m_1p+a_1=m_2q+a_2$，其中 $p, q$ 是整数，则有 $m_1p-m_2q=a_2-a_1$。

由裴蜀定理，当 $a_2-a_1$ 不能被 $\gcd(m_1,m_2)$ 整除时，无解；

其他情况下，可以通过扩展欧几里得算法解出来一组可行解 $(p, q)$；

则原来的两方程组成的模方程组的解为 $x\equiv b\pmod M$，其中 $b=m_1p+a_1$，$M=\text{lcm}(m_1, m_2)$。

多组方程时按以上方法进行合并。

#include<iostream>
#include<algorithm>
typedef __int128 ll;
using namespace std;
ll A, B;
void exgcd(ll a, ll b, ll &d, ll& x, ll& y) {
    if (!b) { d = a;x = 1; y = 0; }
    else { exgcd(b, a % b, d, y, x);y -= x * (a / b); }
}
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;
}
void merge(ll a, ll b) {
    ll p, q, d;
    exgcd(B, b, d, p, q);
    ll c = a - A;
    if (c % d) {
        puts("-1");
        exit(0);
    }
    p = p * (c / d) % (b / d);
    if (p < 0) p += b / d;
    ll mod = lcm(b, B);
    A = (B * p + A) % mod;
    if (A < 0) A += mod;
    B = mod;
}
int main() {
    int t;
    cin >> t;
    while (t--) {
        long long a, b;
        cin >> b >> a;
        if (A == 0 && B == 0)
            A = a, B = b;
        else
            merge(a, b);
    }
    cout << (long long)(A%B) << endl;
    return 0;
}

扩展BSGS

利用扩展BSGS算法解决离散对数问题：$a^x \equiv b \pmod p$

若$a、p$互质，则:

令 $x = A \left \lceil \sqrt p \right \rceil - B$，其中 $0\le A,B \le \left \lceil \sqrt p \right \rceil$，则有 $a^{A\left \lceil \sqrt p \right \rceil -B} \equiv b \pmod p$，稍加变换，则有 $a^{A\left \lceil \sqrt p \right \rceil} \equiv ba^B \pmod p$。

我们已知的是 $a,b$，所以我们可以先算出等式右边的 $ba^B$ 的所有取值，枚举 $B$，用 hash/map 存下来，然后逐一计算 $a^{A\left \lceil \sqrt p \right \rceil}$，枚举 $A$，寻找是否有与之相等的 $ba^B$，从而我们可以得到所有的 $x$，$x=A \left \lceil \sqrt p \right \rceil - B$。

注意到 $A,B$ 均小于 $\left \lceil \sqrt p \right \rceil$，所以时间复杂度为 $\Theta\left (\sqrt p\right )$，用 map 则多一个 $\log$。

若$a、p$不互质，则：

当 $a\perp p$ 时，在模 $p$ 意义下 $a$ 存在逆元，因此可以使用 BSGS 算法求解。于是我们想办法让他们变得互质。

具体地，设 $d_1=\gcd(a,p)$。如果 $d_1\nmid b$，则原方程无解。否则我们把方程同时除以 $d_1$，得到

$$\frac{a}{d_1}\cdot a^{x-1}\equiv \frac{b}{d_1}\pmod{\frac{p}{d_1}}$$

如果 $a$ 和 $\frac{p}{d_1}$ 仍不互质就再除，设 $d_2=\gcd\left(a,\frac{p}{d_1}\right)$。如果 $d_2\nmid \frac{b}{d_1}$，则方程无解；否则同时除以 $d_2$ 得到

$$\frac{a^2}{d_1d_2}\cdot a^{x-2}≡\frac{b}{d_1d_2} \pmod{\frac{p}{d_1d_2}}$$

同理，这样不停的判断下去。直到 $a\perp \frac{p}{d_1d_2\cdots d_k}$。

记 $D=\prod_{i=1}^kd_i$，于是方程就变成了这样：

$$\frac{a^k}{D}\cdot a^{x-k}\equiv\frac{b}{D} \pmod{\frac{p}{D}}$$

由于 $a\perp\frac{p}{D}$，于是推出 $\frac{a^k}{D}\perp \frac{p}{D}$。这样 $\frac{a^k}{D}$ 就有逆元了，于是把它丢到方程右边，这就是一个普通的 BSGS 问题了，于是求解 $x-k$ 后再加上 $k$ 就是原方程的解啦。

注意，不排除解小于等于 $k$ 的情况，所以在消因子之前做一下 $\Theta(k)$ 枚举，直接验证 $a^i\equiv b \pmod p$，这样就能避免这种情况。

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

int BSGS(int a, int b, int p) {
    if (1 % p == b % p) return 0;
    int k = sqrt(p) + 1;
    map<int, int> hs;
    for (int y = 0, r = b % p; y < k; y++) {
        hs[r] = y;
        r = (LL)r * a % p;
    }
    int ak = 1;
    for (int i = 1; i <= k; i++) ak = (LL)ak * a % p;
    for (int x = 1, l = ak; x <= k; x++) {
        if (hs.count(l)) return k * x - hs[l];
        l = (LL)l * ak % p;
    }
    return -INF;
}

int exBSGS(int a, int b, int p) {
    b = (b % p + p) % p;
    if (1 % p == b % p) return 0;
    int x, y;
    int d = exgcd(a, p, x, y);
    if (d > 1) {
        if (b % d) return -INF;
        exgcd(a / d, p / d, x, y);
        return exBSGS(a, (LL)b / d * x % (p / d), p / d) + 1;
    }
    return BSGS(a, b, p);
}