DDOSvoid's Blog

Luogu P4449 于神之怒加强版

发表于 更新于 分类于
本文字数： 1.9k 阅读时长 ≈ 2 分钟

题目描述

https://www.luogu.com.cn/problem/P4449

简要题意：求 $\sum_{i=1}^n\sum_{j=1}^m(i,j)^k$

$n,m\le 5\times 10^6,T\le 2\times 10^3$

Solution

令 $f(n)=\sum_{d|n}d^k\mu(\frac{n}{d})$，容易发现这是一个积性函数，看起来好像不是很能直接线性筛（其实是可以的

然后这东西在 $n=p^c$ 的时候可以简单的求解，所以先算出来 $n=p^c$ 的值，然后再求剩下的值

时间复杂度同样是 $O(n)$，单次询问可以做到 $O(\sqrt n)$

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
#include <iostream>
#include <cstdio>
#define maxn 5000010
#define ll long long 
using namespace std;

const int p = 1000000007;

ll pow_mod(ll x, ll n) {
    ll s = 1;
    for (; n; n >>= 1) {
        if (n & 1) s = s * x % p;
        x = x * x % p;
    }
    return s;
}

int n, m, k; 

int pri[maxn], cnt, a[maxn]; bool isp[maxn]; ll f[maxn]; 
void init_isp(int n) {
    for (int i = 2; i <= n; ++i) {
        if (!isp[i]) pri[++cnt] = i, a[i] = 1;
        for (int j = 1; j <= cnt && i * pri[j] <= n; ++j) {
            isp[i * pri[j]] = 1;
            if (i % pri[j] == 0) { a[i * pri[j]] = a[i]; break; }
            a[i * pri[j]] = i;
        }
    } f[1] = 1; 
    for (int i = 1; i <= cnt; ++i) {
        ll pk = pow_mod(pri[i], k), inv = pow_mod(pk, p - 2); 
        for (ll j = pri[i], t = pk; j <= n; j *= pri[i], t = t * pk % p)
            f[j] = t * (1 - inv) % p; 
    }
    for (int i = 2; i <= n; ++i)
        if (a[i] > 1) f[i] = f[i / a[i]] * f[a[i]] % p; 
    for (int i = 1; i <= n; ++i) f[i] = (f[i] + f[i - 1]) % p; 
}

void work() {
    cin >> n >> m; ll ans = 0; if (n > m) swap(n, m); 
    for (int l = 1, r; l <= n; l = r + 1) {
        r = min(n / (n / l), m / (m / l));
        ans = (ans + (f[r] - f[l - 1]) * (n / l) % p * (m / l)) % p;
    } cout << (ans + p) % p << "\n"; 
} 

int main() {
    ios::sync_with_stdio(false);
    cin.tie(nullptr); cout.tie(nullptr);

    int T; cin >> T >> k; init_isp(5000000);
    while (T--) work(); 
    return 0;
}