BZOJ4407 于神之怒加强版 <莫比乌斯反演>

Problem

于神之怒加强版

$\mathrm{Time\;Limit:\;80\;Sec}$
$\mathrm{Memory\;Limit:\;512\;MB}$

Description

给下$N,M,K$，计算$\sum_{i=1}^{n}\sum_{j=1}^{m}\gcd(i,j)^k\mod(10^9+7)$的值。

Input

输入有多组数据，输入数据的第一行两个正整数$T,K$，代表有$T$组数据，$K$的意义如上所示，下面第$2$行到第$T+1$行，每行为两个正整数$N,M$，其意义如上式所示。

Output

对于每一个询问，输出一行一个数作为回答。

Sample Input

1 2
3 3

Sample Output

20

HINT

$1\le N,M,K\le5\times 10^6,\;1\le T\le2000$
官方题解

标签：莫比乌斯反演

Solution

先套路转换出$\mu$：

$$\begin{aligned} Ans&=\sum_{i=1}^{n}\sum_{j=1}^{m}\gcd(i,j)^k\\ &=\sum_{d=1}^{\min(n,m)}d^k\sum_{i=1}^{n}\sum_{j=1}^{m}[\gcd(i,j)=d]\\ &=\sum_{d=1}^{\min(n,m)}d^k\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{m}{d}\rfloor}\sum_{x|\gcd(i ,j)}\mu(x)\\ &=\sum_{d=1}^{\min(n,m)}d^k\sum_{x=1}^{\lfloor\frac{\min(n,m)}{d}\rfloor}\mu(x)\lfloor\frac{n}{x\cdot d}\rfloor\lfloor\frac{m}{x\cdot d}\rfloor\\ &=\sum_{t=1}^{\min(n,m)}\sum_{d=1}^{t}d^k\mu(\frac{t}{d})\lfloor\frac{n}{t}\rfloor\lfloor\frac{m}{t}\rfloor\\ &=\sum_{t=1}^{\min(n,m)}\lfloor\frac{n}{t}\rfloor\lfloor\frac{m}{t}\rfloor\sum_{d=1}^{t}d^k\mu(\frac{t}{d})\\ \end{aligned}$$

那么每次询问对于前半部分可以根号分块，随后需要$O(1)$计算后半部分的值，因而需要线筛预处理后半部分的值。

令$f(t)=\sum_{d=1}^{t}d^k\mu(\frac{t}{d})$，$t=p_1^{a_1}\times p_2^{q_2}\times\cdots p_l^{a_l}$
由于$y=d^k$和$y=\mu(d)$均为积性函数，因而$f$也为积性函数。
那么就有

$$\begin{aligned} f(t)&=\prod_{i=1}^{l}f(p_i^{a_i})\\&=\prod_{i=1}^{l}\sum_{j=0}^{a_i}p_i^{j\cdot k}\mu(\frac{p_i^{a_i}}{p_i^j})\\ &=\prod_{i=1}^{l}\sum_{j=0}^{a_i}p_i^{j\cdot k}\mu(p_i^{a_i-j})\\ \end{aligned}$$

易知当$j\in[0,a_i-1)$时，均有$\mu(p_i^{a_i-j})=0$，因此有

$$\begin{aligned} f(t)&=\prod_{i=1}^{l}\sum_{j=0}^{a_i}p_i^{j\cdot k}\mu(p_i^{a_i-j})\\ &=\prod_{i=1}^{l}(\mu(p_i)\cdot p_i^{k\cdot(a_i-1)}+\mu(1)\cdot p_i^{k\cdot a_i})\\ &=\prod_{i=1}^{l}(p_i^{k\cdot a_i}-p_i^{k\cdot(a_i-1)})\\ &=\prod_{i=1}^{l}p_i^{k\cdot(a_i-1)}(p_i-1) \end{aligned}$$

接下来考虑在线筛中如何处理。
首先，对于所有质数，均有$f(p)=p^k-1$。
而对于合数，假设当前筛到的数是$x$，对于一个比它小的素数$p$，有两种情况：

1. 若$p|x$，设$p$在$x$分解质因数中的次数为$a$，那么$f(x\times p)$相比于$f(x)$而言，在含$p$的约数中$p$的次数都增加了$1$，否则$\mu=0$无贡献。因此$k\cdot(a-1)$增加了$k$，这样总共扩大了$p^k$倍，故$f(x\times p)=f(x)\times p^k$；
2. 若$p\nmid x$，由积性函数可知$f(x\times p)=f(x)\times f(p)$。

这样就可以$O(n\log{n})$筛出$f$的函数值，每次询问$O(T\sqrt{n})$根号分块，总复杂度$O(n\log{n}+T\sqrt{n})$。

Code

#include <bits/stdc++.h>
#define MAX_N 5000000
#define MOD 1000000007
using namespace std;
typedef long long lnt;
template <class T> inline void read(T &x) {
    x = 0; int c = getchar(), f = 1;
    for (; !isdigit(c); c = getchar()) if (c == 45) f = -1;
    for (; isdigit(c); c = getchar()) (x *= 10) += f*(c-'0');
}
lnt k, cnt, ans, f[MAX_N+5], p[MAX_N+5], pri[MAX_N+5];    bool NotPri[MAX_N+5];
lnt PM(lnt x, lnt y) {if (!y) return 1LL; lnt ret = PM(x, y>>1); return (y&1) ? ret*ret%MOD*x%MOD : ret*ret%MOD;}
void init() {
    NotPri[1] = true, f[1] = 1;
    for (int i = 2; i <= MAX_N; i++) {
        if (!NotPri[i]) pri[cnt++] = i, p[i] = PM(i, k), f[i] = p[i]-1;
        for (int j = 0; j < cnt; j++) {
            if (i*pri[j] > MAX_N) break; NotPri[i*pri[j]] = true;
            if (i%pri[j]) f[i*pri[j]] = f[i]*f[pri[j]]%MOD;
            else {f[i*pri[j]] = f[i]*p[pri[j]]%MOD; break;}
        }
    }
    for (int i = 2; i <= MAX_N; i++) (f[i] += f[i-1]) %= MOD;
}
int main() {
    int T;	read(T), read(k), init();
    while (T--) {
        lnt n, m;	read(n), read(m), ans = 0;
        for (lnt l = 1, r; l <= min(n, m); l = r+1)
            r = min(n/(n/l), m/(m/l)), (ans += (n/l)*(m/l)%MOD*(f[r]-f[l-1]+MOD)%MOD) %= MOD;
        printf("%lld\n", ans);
    }
    return 0;
}