bzoj4407 于神之怒加强版

 传送门：http://www.lydsy.com/JudgeOnline/problem.php?id=4407

【题解】

推一波公式

你还是需要前置技能：

那么好像可以O(Tnlogn)直接暴力啊！

当然是两遍根号分块变成O(Tn)啊

好消息过不去。

当你莫比乌斯反演发现复杂度不对的话怎么办？继续瞎**化简！

令x=pd，那么有

哎这个式子看起来很和善

后面那个好像是积性函数（逃

那么线性筛就行啦！！！

问题是怎么筛呢。我们把每个因数分开考虑后发现。

如果这个因数i只出现一次，贡献显然为（i的k次方-1），这可以直接写出来的qwq

如果出现了多次，我们发现只有mu(1)和mu(i)是有值的，所以只考虑这两项，写出来发现是在前面的基础上乘上了（i的k次方）

于是就能筛啦！

# include <stdio.h>
# include <string.h>
# include <algorithm>
// # include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef unsigned long long ull;
const int M = 5e5 + 10;
const int mod = 1e9+7;

# define RG register
# define ST static

int n, m, K;

const int F = 5000000;
bool isnp[F + 10];
int p[F/3], pn=0;
int f[F + 10], s[F/3];

inline int pwr(int a, int b) {
    int ret = 1;
    while(b) {
        if(b&1) ret = 1ll * ret * a % mod;
        a = 1ll * a * a % mod;
        b >>= 1;
    }
    return ret;
}

inline void sieve() {
    f[1] = 1;  
    for (int i=2; i<=F; ++i) {
        if(!isnp[i]) {
            p[++pn] = i;
            s[pn] = pwr(i, K); 
            f[i] = s[pn] - 1; 
            if(f[i] < 0) f[i] += mod; 
        }
        for (int j=1; j<=pn && i*p[j]<=F; ++j) {
            isnp[i*p[j]] = 1;
            if(i%p[j] == 0) {
                f[i*p[j]] = 1ll * f[i] * s[j] % mod;
                break;
            }
            f[i*p[j]] = 1ll * f[i] * f[p[j]] % mod;
        }
    }
    for (int i=1; i<=F; ++i) {
        f[i] = f[i] + f[i-1];
        if(f[i] >= mod) f[i] -= mod;
    }
}

inline void sol() {
    scanf("%d%d", &n, &m);    
    if(n>m) swap(n, m);
    int ans = 0;
    for (int i=1, lst; i<=n; i=lst+1) {
        lst = min(n/(n/i), m/(m/i));
        ans = ans + 1ll * (f[lst] + mod - f[i-1]) * (n/i) % mod * (m/i) % mod;
        if(ans >= mod) ans -= mod;
    } 
    printf("%d
", ans);     
}

int main() {
    int T; scanf("%d%d", &T, &K); 
    sieve(); 
    while(T--) sol(); 
    return 0;
}