BZOJ 4407 于神之怒加强版

【题目描述】#

给定n,m,k,计算 $\sum_{i=1}^{n}{\sum_{j=1}^{m}{gcd(i,j)^k}}$ 对1000000007取模的结果

【输入格式】#

多组数据。第一行是两个数T,K; 之后的T行，每行两个整数n，m；

【输出格式】#

K行，每行一个结果

【样例输入】#

1 2
3 3

【样例输出】#

20

【提示】#

T<=2000,1<=N,M,K<=5000000

题解#

第一步推式子

$\sum_{i=1}^{n}{\sum_{j=1}^{m}{gcd(i,j)^k}}$
$=\sum_{d=1}^{n}{d^k\sum_{i=1}^{n}{\sum_{j=1}^{m}{[gcd(i,j)==d]}}}$
$=\sum_{d=1}^{n}{d^k\sum_{i=1}^{\lfloor \frac{n}{d} \rfloor}{\sum_{j=1}^{\lfloor \frac{m}{d} \rfloor}{[gcd(i,j)==1]}}}$
$=\sum_{d=1}^{n}{d^k\sum_{i=1}^{\lfloor \frac{n}{d} \rfloor}{\sum_{j=1}^{\lfloor \frac{m}{d} \rfloor}{\sum_{c|i,c|j}{\mu{c}}}}}$
另 $T = dc$
$=\sum_{T=1}^{n}{\sum_{d|T}{d^k \mu{\frac{T}{d}} \lfloor \frac{n}{T} \rfloor \lfloor \frac{m}{T} \rfloor}}$
$=\sum_{T=1}^{n}{\lfloor \frac{n}{T} \rfloor \lfloor \frac{m}{T} \rfloor \sum_{d|T}{d^k \mu{\frac{T}{d}}}}$

然后我们需要线性筛出 $\sum_{d|T}{d^k \mu{\frac{T}{d}}}$ 就可以了
设 $g(T) = \sum_{d|T}{d^k \mu{\frac{T}{d}}}$
显然是积性函数别问我怎么知道的
当互质时 $g(T*P) = g(T) * g(p)$
当不互质时又要推狮子
根据 $\mu$ 函数的定义
当且仅当不含平方因子时不为零
所以质因子只有选和不选两种状态，有用的状态数是不变的
每一个有用的 $\mu$ 都乘了一个p^k; 所以 $g(T*p) = g(T) * p^k$

1
#define LL long long
2

3
#include <cstdio>
4
#include <cstring>
5
#include <algorithm>
6
using namespace std;
7
int prime[500005], tot;
8
const int MOD = 1000000007;
9
bool isnprime[5000005];
10
LL g[5000005];
11
LL primeK[500005], n, k, m;
12
LL N = 5000000;
13
LL pow_mod(LL a, LL b)
14
{
15
    LL ans = 1;
16
    while (b)
17
    {
18
        if (b & 1)
19
            ans = ans * a % MOD;
20
        b >>= 1;
21
        a = a * a % MOD;
22
    }
23
    return ans;
24
}
25
void Get_g()
26
{
27
    g[1] = 1;
28
    for (int i = 2; i <= N; i++)
29
    {
30
        if (!isnprime[i])
31
        {
32
            prime[++tot] = i;
33
            g[i] = (primeK[tot] = pow_mod(i, k)) - 1;
34
        }
35
        for (int j = 1; j <= tot; j++)
36
        {
37
            if (i * prime[j] > N)
38
                break;
39
            isnprime[i * prime[j]] = 1;
40
            if (i % prime[j] == 0)
41
            {
42
                g[i * prime[j]] = g[i] * primeK[j] % MOD;
43
                break;
44
            }
45
            g[i * prime[j]] = g[i] * g[prime[j]] % MOD;
46
        }
47
    }
48
    for (int i = 2; i <= N; i++)
49
        g[i] += g[i - 1], g[i] %= MOD;
50
}
51
int main()
52
{
53
    freopen("bzoj_4407.in","r",stdin);
54
    freopen("bzoj_4407.out","w",stdout);
55
    int t;
56
    scanf("%d%lld", &t, &k);
57
    Get_g();
58
    while (t--)
59
    {
60
        scanf("%lld%lld", &n, &m);
61
        int last;
62
        LL ret = 0;
63
        if (n > m)
64
            swap(n, m);
65
        for (int i = 1; i <= n; i = last + 1)
66
        {
67
            last = min(n / (n / i), m / (m / i));
68
            (ret += (n / i) * (m / i) % MOD * (g[last] - g[i - 1] + MOD) % MOD)%=MOD;
69
        }
70
        printf("%lld\n", ret);
71
    }
72
}