BZOJ4407 于神之怒加强版 - 莫比乌斯反演
题解
非常裸的莫比乌斯反演。
但是反演完还需要快速计算一个积性函数（我直接用$nlogn$卷积被TLE了
我也不想再写一遍了
代码
1 #include<cstring> 2 #include<cstdio> 3 #include<algorithm> 4 #define ll long long 5 #define rd read() 6 using namespace std; 7 8 const ll mod = 1e9 + 7; 9 10 const int N = 5e6; 11 12 int T, n, m, k; 13 int mu[N], pri[N], tot, vis[N], low[N]; 14 ll sum[N]; 15 16 int read() { 17 int X = 0, p = 1; char c = getchar(); 18 for(; c > '9' || c < '0'; c = getchar()) if(c == '-') p = -1; 19 for(; c >= '0' && c <= '9'; c = getchar()) X = X * 10 + c - '0'; 20 return X * p; 21 } 22 23 ll fpow(ll a, ll b) { 24 ll re = 1; 25 for(; b; b >>= 1, a = a * a % mod) if(b & 1) re = re * a % mod; 26 return re; 27 } 28 29 void init() { 30 mu[1] = 1; 31 sum[1] = 1; 32 for(int i = 2; i < N; ++i) { 33 if(!vis[i]) mu[i] = -1, pri[++tot] = i, sum[i] = fpow(i, k) - 1, low[i] = i; 34 for(int j = 1; j <= tot && i * pri[j] < N; ++j) { 35 vis[i * pri[j]] = 1; 36 if(i % pri[j] == 0) { 37 low[i * pri[j]] = low[i] * pri[j]; 38 mu[i * pri[j]] = 0; 39 if(low[i] == i) sum[i * pri[j]] = sum[i] * fpow(pri[j], k) % mod; 40 else sum[i * pri[j]] = sum[i / low[i]] * sum[pri[j] * low[i]] % mod; 41 break; 42 } 43 low[i * pri[j]] = pri[j]; 44 mu[i * pri[j]] = -mu[i]; 45 sum[i * pri[j]] = sum[i] * sum[pri[j]] % mod; 46 } 47 } 48 //for(ll i = 1; i < N; ++i) { 49 // ll tmp = fpow(i, k); 50 // for(ll j = 1; j < N && i * j < N; ++j) sum[i * j] = (sum[i * j] + tmp * mu[j]) % mod; 51 //} 52 for(int i = 1; i < N; ++i) sum[i] = (sum[i - 1] + sum[i]) % mod; 53 } 54 55 ll cal() { 56 ll ans = 0; 57 for(int i = 1, j = 1; i <= min(n, m); i = j + 1) { 58 j = n / (n / i); 59 j = min(j, m / (m / i)); 60 ans += (sum[j] - sum[i - 1]) % mod * (n / i) % mod * (m / i) % mod; 61 ans %= mod; 62 } 63 ans = (ans % mod + mod) % mod; 64 return ans; 65 } 66 67 int main() 68 { 69 T = rd; k = rd; 70 init(); 71 for(; T; T--) { 72 n = rd; m = rd; 73 printf("%lld\n", cal()); 74 } 75 }
