【题解】P4449 于神之怒加强版

不太会线筛求函数前缀和，记一下。

思路

莫反 + 线筛。

首先根据莫反套路有：

\[\begin{aligned} &\sum\limits_{i = 1}^n \sum\limits_{j = 1}^m \gcd(i, j)^k \\ &= \sum\limits_{d = 1}^{\min(n, m)} d^k \sum\limits_{i = 1}^{\lfloor \frac{n}{d} \rfloor} \sum\limits_{j = 1}^{\lfloor \frac{m}{d} \rfloor} [\gcd(i, j) = 1] \\ &= \sum\limits_{d = 1}^{\min(n, m)} d^k \sum\limits_{i = 1}^{\lfloor \frac{n}{d} \rfloor} \sum\limits_{j = 1}^{\lfloor \frac{m}{d} \rfloor} \sum\limits_{x \mid \gcd(i, j)} \mu(x) \\ &= \sum\limits_{d = 1}^{\min(n, m)} d^k \sum\limits_{x = 1}^{\min(\lfloor \frac{n}{d} \rfloor, \lfloor \frac{m}{d} \rfloor)} \mu(x) \cdot \lfloor \frac{n}{dx} \rfloor \lfloor \frac{m}{dx} \rfloor \\ &= \sum\limits_{T = 1}^{\min(n, m)} \lfloor \frac{n}{T} \rfloor \lfloor \frac{m}{T} \rfloor \sum\limits_{d \mid T} d^k \mu(\frac{T}{d}) \end{aligned} \]

前半部分只需用数论分块处理，现在需要考虑求后半部分。

如果不用数据结构的话，只能用线筛维护后半部分。

令 \(f(x) = \sum\limits_{d \mid x} d^k \mu(\frac{x}{d})\)，考虑线筛求 \(f\).

先考虑 \(f\) 在某一质数 \(p\) 的整次幂 \(p^c\) 处的取值。

根据定义有 \(f(p^c) = \sum\limits_{d \mid p^c} d^k \cdot \mu(\frac{p^c}{d})\)

因为对于含平方因子的数 \(x\)，\(\mu(x) = 0\). 所以 \(f(p^c)\) 可以写成 \(f(p^c) = (p^{c - 1})^k \cdot \mu(p) + (p^c)^k \cdot \mu(1) = p^{ck} - p^{(c - 1)k} = p^{(c - 1)k} \cdot (p^k - 1)\)

所以有 \(f(p^{c + 1}) = p^{ck} \cdot (p^k - 1) = p^k \cdot p^{(c - 1)k} \cdot (p^k - 1) = f(p^c) \cdot p^k\)

最后考虑边界条件 \(f(p)\).

显然有 \(f(p) = \sum\limits_{d \mid p} d^k \mu(\frac{p}{d}) = 1^k \cdot \mu(p) + p^{k} \cdot \mu(1) = p^k - 1\)

所以可以推出线筛的流程：

假设当前已知 \(f(x), f(y)\)，其中 \(x \in \mathbb{N^+}, y \in \mathbb{P}\)，考虑筛出 \(f(xy)\).

1. \(y \mid x\)，令 \(x\) 中 \(y\) 的幂次为 \(c\). 根据积性函数的定义有 \(f(xy) = f(\frac{xy}{y^{c + 1}}) \cdot f(y^{c + 1})\)

   又因为 \(f(y^{c + 1}) = f(y^c) \cdot y^k\)

   所以有 \(f(xy) = f(\frac{x}{y^c}) \cdot f(y^c) \cdot y^k = f(\frac{x}{y^c} \cdot y^c) \cdot y^k = f(x) \cdot y^k\)

2. \(y \nmid x\)，由积性直接可得 \(f(xy) = f(x) \cdot f(y)\)

所以直接线性筛再套一层数论分块即可。

时间复杂度 \(O(T \sqrt{n} + n \log V)\)

另外 \(k = 1\) 时 \(f\) 等价于欧拉函数 \(\varphi\)

代码

#include <cstdio>
#include <vector>
#include <bitset>
using namespace std;

const int maxn = 5e6 + 5;
const int mod = 1e9 + 7;

int t, n, m, k;
int f[maxn];
vector<int> p;
bitset<maxn> vis;

int qpow(int base, int power)
{
    int res = 1;
    while (power)
    {
        if (power & 1) res = 1ll * res * base % mod;
        base = 1ll * base * base % mod;
        power >>= 1;
    }
    return res;
}

void init()
{
    vis[1] = true, f[1] = 1;
    for (int i = 2; i <= 5e6; i++)
    {
        if (!vis[i]) p.push_back(i), f[i] = (qpow(i, k) - 1 + mod) % mod;
        for (int j = 0; (j < p.size()) && (i * p[j] <= 5e6); j++)
        {
            vis[i * p[j]] = true;
            if (i % p[j] == 0) { f[i * p[j]] = 1ll * f[i] * qpow(p[j], k) % mod; break; }
            else f[i * p[j]] = 1ll * f[i] * f[p[j]] % mod;
        }
    }
    for (int i = 2; i <= 5e6; i++) f[i] = (f[i] + f[i - 1]) % mod;
}

int main()
{
    scanf("%d%d", &t, &k);
    init();
    while (t--)
    {
        scanf("%d%d", &n, &m);
        if (n > m) swap(n, m);
        int ans = 0;
        for (int l = 1, r; l <= n; l = r + 1)
        {
            r = min(n / (n / l), m / (m / l));
            ans = (ans + 1ll * (f[r] - f[l - 1] + mod) % mod * (n / l) % mod * (m / l) % mod) % mod;
        }
        printf("%d\n", ans);
    }
    return 0;
}