莫比乌斯反演

  1. $F(n)=\sum_{d|n}f(d) \Rightarrow f(n)=\sum_{d|n}\mu(d)F(\frac{n}{d})$

  2. $F(n)=\sum_{n|d}f(d) \Rightarrow f(n)=\sum_{n|d}\mu(\frac{d}{n})F(d)$

[POI2007]Zap

题目链接：bzoj - 1101

题意：对于给定的整数a，b和d，有多少正整数对x，y，满足x$\leq$a，y$\leq$b，并且gcd(x,y)=d

思路：设$F(n)$表示gcd(x,y)为n的倍数(包括n本身)的数量，$f(n)$表示gcd(x,y)为n的数量

显然有$F(n)=\sum_{n|d}f(d)$

所以由莫比乌斯反演得$f(n)=\sum_{n|d}\mu(\frac{d}{n})F(d)$，而$F(d)=\left \lfloor \frac{a}{d} \right \rfloor \left \lfloor \frac{b}{d} \right \rfloor$(根据$F(n)$的定义理解)

代入$f(n)=\sum_{n|d}\mu(\frac{d}{n}) \left \lfloor \frac{a}{d} \right \rfloor \left \lfloor \frac{b}{d} \right \rfloor$

令$k=\frac{d}{n}$可得$f(n)=\sum_{k=1}^{min(\left \lfloor \frac{a}{n} \right \rfloor , \left \lfloor \frac{b}{n} \right \rfloor)}\mu(k) \left \lfloor \frac{a}{kn} \right \rfloor \left \lfloor \frac{b}{kn} \right \rfloor$

预处理出$\mu(k)$的前缀和，然后整除分块计算即可

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
#include <iostream>
#include <algorithm>
#include <cstring>
#include <cstdio>

using namespace std;

const int N = 500010;

typedef long long ll;

int isp[N], p[N], mu[N], sum[N];
int tot, T, a, b, d;

void init(int n)
{
    mu[1] = 1;
    for (int i = 2; i <= n; i++) {
        if (!isp[i]) {
            p[++tot] = i;
            mu[i] = -1;
        }
        for (int j = 1; j <= tot && p[j] <= n / i; j++) {
            isp[i * p[j]] = 1;
            if (0 == i % p[j]) {
                mu[i * p[j]] = 0;
                break;   
            }
            else mu[i * p[j]] = -mu[i];
        }
    }
    for (int i = 1; i <= n; i++) sum[i] = sum[i - 1] + mu[i];
}

int main()
{
    init(N - 5);
    scanf("%d", &T);
    while (T--) {
        scanf("%d%d%d", &a, &b, &d);
        ll res = 0;
        for (int l = 1, r = 1; l <= min(a / d, b / d); l = r + 1) {
            r = min(a / (a / l), b / (b / l));
            res += 1ll * (sum[r] - sum[l - 1]) * (a / (l * d)) * (b / (l * d));
        }
        printf("%lld\n", res);
    }
    return 0;
}

YY的GCD

题目链接：bzoj - 2820

题意：对于给定的整数n，m，有多少正整数对x，y，满足x$\leq$n，y$\leq$m，并且gcd(x,y)为质数

思路：设$F(n)$表示gcd(x,y)为n的倍数(包括n本身)的数量，$f(n)$表示gcd(x,y)为n的数量

显然有$F(n)=\sum_{n|d}f(d)$

所以由莫比乌斯反演得$f(n)=\sum_{n|d}\mu(\frac{d}{n})F(d)$，而$F(d)=\left \lfloor \frac{a}{d} \right \rfloor \left \lfloor \frac{b}{d} \right \rfloor$(根据$F(n)$的定义理解)

代入$f(n)=\sum_{n|d}\mu(\frac{d}{n}) \left \lfloor \frac{a}{d} \right \rfloor \left \lfloor \frac{b}{d} \right \rfloor$

所以$res=\sum_{p\in prime}f(p)=\sum_{p\in prime}\sum_{p|d}\mu(\frac{d}{p}) \left \lfloor \frac{a}{d} \right \rfloor \left \lfloor \frac{b}{d} \right \rfloor$

因为对于质数p的每个倍数d求$\mu(\frac{d}{p}) \left \lfloor \frac{a}{d} \right \rfloor \left \lfloor \frac{b}{d} \right \rfloor$所得到的值和对于d的每个质因子p求$\mu(\frac{d}{p}) \left \lfloor \frac{a}{d} \right \rfloor \left \lfloor \frac{b}{d} \right \rfloor$所得到的值是一样的

所以$res=\sum_{d=1}^{min(n,m)} \sum_{p\in{prime}\ p|d} \mu(\frac{d}{p}) \left \lfloor \frac{a}{d} \right \rfloor \left \lfloor \frac{b}{d} \right \rfloor$

用整数分块化简后得$\sum_{d=1}^{min(n,m)}(\left \lfloor \frac{a}{d} \right \rfloor \left \lfloor \frac{b}{d} \right \rfloor \sum_{p\in{prime}\ p|d} \mu(\frac{d}{p}))$

令$c(d)=\sum_{p\in{prime}\ p|d} \mu(\frac{d}{p})$

预处理出$c(d)$的前缀和$sum(d)$，利用整数分块计算即可

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
#include <iostream>
#include <algorithm>
#include <cstring>
#include <cstdio>

using namespace std;

typedef long long ll;

const int N = 10000010;

int isp[N], p[N], mu[N];
int T, tot, n, m;
ll sum[N];

void init(int n)
{
    mu[1] = 1;
    for (int i = 2; i <= n; i++) {
        if (!isp[i]) {
            p[++tot] = i;
            mu[i] = -1;
        }
        for (int j = 1; j <= tot && p[j] <= n / i; j++) {
            isp[i * p[j]] = 1;
            if (0 == i % p[j]) {
                mu[i * p[j]] = 0;
                break;
            }
            else mu[i * p[j]] = -mu[i];
        }
    }
    for (int i = 1; i <= tot; i++)
        for (int j = 1; p[i] * j <= n; j++)
            sum[p[i] * j] += mu[j];
    for (int i = 1; i <= n; i++) sum[i] += sum[i - 1];
}

int main()
{
    init(N - 5);
    scanf("%d", &T);
    while (T--) {
        scanf("%d%d", &n, &m);
        ll res = 0;
        for (int l = 1, r = 0; l <= min(n, m); l = r + 1) {
            r = min(n / (n / l), m / (m / l));
            res = res + 1ll * (sum[r] - sum[l - 1]) * (n / l) * (m / l);
        }
        printf("%lld\n", res);
    }
    return 0;
}

[HAOI2011]Problem b

题目链接：bzoj - 2301

题意：对于给定的整数a，b，c，d和k，有多少正整数对x，y，满足a$\leq$x$\leq$b，c$\leq$y$\leq$d，并且gcd(x,y)=k

思路：利用第1题的结果和容斥原理即可得到答案，$res=solve(b,d)-solve(a-1,d)-solve(b,c-1)+solve(a-1,c-1)$

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
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
#include <iostream>
#include <algorithm>
#include <cstring>
#include <cstdio>

using namespace std;

const int N = 50010;

int isp[N], p[N], mu[N], sum[N];
int tot, T, a, b, c, d, k;

inline int read()
{
    int s = 0, f = 1;
    char ch = getchar();
    while (ch < '0' || ch > '9') {
        if ('-' == ch) f = -1;
        ch = getchar();
    }
    while (ch >= '0' && ch <= '9') {
        s = s * 10 + ch - '0';
        ch = getchar();
    }
    return s * f;
}

inline void write(int x)
{
    if (x < 0) {
        putchar('-');
        x = -x;
    }
    if (x > 9) write(x / 10);
    putchar(x % 10 + '0');
}

void init(int n)
{
    mu[1] = 1;
    for (int i = 2; i <= n; i++) {
        if (!isp[i]) {
            p[++tot] = i;
            mu[i] = -1;
        }
        for (int j = 1; j <= tot && p[j] <= n / i; j++) {
            isp[i * p[j]] = 1;
            if (0 == i % p[j]) {
                mu[i * p[j]] = 0;
                break;   
            }
            else mu[i * p[j]] = -mu[i];
        }
    }
    for (int i = 1; i <= n; i++) sum[i] = sum[i - 1] + mu[i];
}

int solve(int a, int b, int n)
{
    int res = 0, ta = a / n, tb = b / n;
    for (int l = 1, r = 1; l <= min(ta, tb); l = r + 1) {
        r = min(a / (a / l), b / (b / l));
        res += (sum[r] - sum[l - 1]) * (a / (n * l)) * (b / (n * l));
    }
    return res;
}

int main()
{
    init(N - 5);
    T = read();
    while (T--) {
        a = read(), b = read(), c = read(), d = read(), k = read();
        int res = solve(b, d, k) - solve(a - 1, d, k) - solve(b, c - 1, k) + solve(a - 1, c - 1, k);
        write(res);
        puts("");
    }
    return 0;
}

[Sdoi2014]数表

题目链接：bzoj - 3529

题意：有一张n*m的数表，其第i行第j列的数值为能同时整除i和j的所有自然数之和，给定a，计算表中不大于a的数之和

思路：假设n<m，先不考虑不大于a的这个约束条件，那么$res=\sum_{i=1}^{n}\sum_{j=1}^{m}\sigma (gcd(i,j))$，其中$\sigma $表示约数之和

令$d=gcd(i,j)$，枚举d可得$res=\sum_{d=1}^{n}\sum_{i=1}^{n}\sum_{j=1}^{m}\sigma(d)[gcd(i,j)=d]$

将$\sigma(d)$提前得$res=\sum_{d=1}^{n}[\sigma(d)\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(i,j)=d]]$

后半部分$\sum_{i=1}^{n}\sum_{j=1}^{m}[gcd(i,j)=d]$相当于求有多少对i，j满足i$\leq n$，$j\leq m$且$gcd(i,j)=d$

后半部分用莫比乌斯化简一下得$res=\sum_{d=1}^{n}[\sigma(d)\sum_{k=1}^{\left \lfloor \frac{n}{d} \right \rfloor} \mu(k) \left \lfloor \frac{n}{kd} \right \rfloor \left \lfloor \frac{m}{kd} \right \rfloor]$

令$T=kd$得$res=\sum_{d=1}^{n}[\sigma(d)\sum_{T=1\ d|T}^{n} \mu(\frac{T}{d}) \left \lfloor \frac{n}{T} \right \rfloor \left \lfloor \frac{m}{T} \right \rfloor]$

交换一下顺序得$res=\sum_{T=1}^{n}[\left \lfloor \frac{n}{T} \right \rfloor \left \lfloor \frac{m}{T} \right \rfloor \sum_{d=1\ d|T}^{n}\sigma(d) \mu(\frac{T}{d})]$

现在考虑限制条件，只有当$\sigma(d)\leq a$时，$\sigma(d) \mu(\frac{T}{d})$才对答案有贡献，所以我们先预处理出$\sigma(d)$，对询问按a排序，对$\sigma(d)$排序，当a变大时，某些$\sigma(d) \mu(\frac{T}{d})$开始对答案产生贡献，此时枚举d得倍数T，然后用树状数组维护$\sigma(d) \mu(\frac{T}{d})$的前缀和就可以计算出$\sum_{d=1\ d|T}^{n}\sigma(d) \mu(\frac{T}{d})$，这个d对别的询问也会产生贡献(询问是按a排序的)，对于每个询问用整数分块计算

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
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
#include <iostream>
#include <algorithm>
#include <cstring>
#include <cstdio>

using namespace std;

const int N = 100010;

struct node {
    int n, m, a, id;
};

int T, c[N], mu[N], res[N];
int isp[N], p[N], tot;
node q[N], d[N];

bool cmp(node a, node b)
{
    return a.a < b.a;
}


void init(int n)
{
    mu[1] = 1;
    for (int i = 2; i <= n; i++) {
        if (!isp[i]) {
            p[++tot] = i;
            mu[i] = -1;
        }
        for (int j = 1; j <= tot && p[j] <= n / i; j++) {
            isp[i * p[j]] = 1;
            if (0 == i % p[j]) {
                mu[i * p[j]] = 0;
                break;
            }
            mu[i * p[j]] = -mu[i];
        }
    }
    for (int i = 1; i <= n; i++) {
        d[i].id = i;
        for (int j = i; j <= n; j += i) {
            d[j].a += i;
        }
    }
    sort(d + 1, d + n + 1, cmp);
}

inline int lowbit(int x)
{
    return x & (-x);
}

void add(int x, int v)
{
    while (x <= N - 5) {
        c[x] += v;
        x += lowbit(x);
    }
}

int ask(int x)
{
    int res = 0;
    while (x > 0) {
        res += c[x];
        x -= lowbit(x);
    }
    return res;
}

int solve(int n, int m)
{
    int res = 0;
    for (int l = 1, r = 1; l <= n; l = r + 1) {
        r = min(n / (n / l), m / (m / l));
        res += (n / l) * (m / l) * (ask(r) - ask(l - 1));
    }
    return res;
}

int main()
{
    init(N - 5);
    scanf("%d", &T);
    for (int i = 1; i <= T; i++) {
        scanf("%d%d%d", &q[i].n, &q[i].m, &q[i].a);
        if (q[i].n > q[i].m) swap(q[i].n, q[i].m);
        q[i].id = i;
    }
    sort(q + 1, q + T + 1, cmp);
    for (int i = 1, j = 1; i <= T; i++) {
        while (j <= N - 5 && d[j].a <= q[i].a) {
            for (int k = d[j].id; k <= N - 5; k += d[j].id) {
                add(k, d[j].a * mu[k / d[j].id]);
            }
            j++;
        }
        res[q[i].id] = solve(q[i].n, q[i].m);
    }
    for (int i = 1; i <= T; i++) printf("%d\n", res[i] & (~(1 << 31)));
    return 0;
}

Crash的数字表格

题目链接：bzoj - 2154

题意：有一张n*m的数表，其第i行第j列的数值为lcm(i,j)，求数表内所有数的和

思路：假设n<m，$res=\sum_{i=1}^{n} \sum_{j=1}^{m} lcm(i,j)=\sum_{i=1}^{n}\sum_{j=1}^{m}\frac{ij}{gcd(i,j)}$

令$d=gcd(i,j)$，枚举d得$res=\sum_{d=1}^{n}\sum_{i=1}^{n} \sum_{j=1}^{m} \frac{ij}{d} [gcd(i,j)=d]=\sum_{d=1}^{n}\sum_{i=1}^{\left \lfloor \frac{n}{d} \right \rfloor} \sum_{j=1}^{\left \lfloor \frac{m}{d} \right \rfloor} ijd [gcd(i,j)=1]$

把d移到前面有$res=\sum_{d=1}^{n}[d\sum_{i=1}^{\left \lfloor \frac{n}{d} \right \rfloor} \sum_{j=1}^{\left \lfloor \frac{m}{d} \right \rfloor} ij [gcd(i,j)=1]]$

根据莫比乌斯得性质得$res=\sum_{d=1}^{n}[d\sum_{i=1}^{\left \lfloor \frac{n}{d} \right \rfloor} \sum_{j=1}^{\left \lfloor \frac{m}{d} \right \rfloor} ij \sum_{k|gcd(i,j)}\mu(k)]$

枚举k得$res=\sum_{d=1}^{n}[d\sum_{k=1}^{\left \lfloor \frac{n}{d} \right \rfloor} \sum_{i=1}^{\left \lfloor \frac{n}{dk} \right \rfloor} \sum_{j=1}^{\left \lfloor \frac{m}{dk} \right \rfloor} ijk^{2}\mu(k)]$

移一下项得$res=\sum_{d=1}^{n}[d\sum_{k=1}^{\left \lfloor \frac{n}{d} \right \rfloor} [ k^{2}\mu(k) \sum_{i=1}^{\left \lfloor \frac{n}{dk} \right \rfloor} i \sum_{j=1}^{\left \lfloor \frac{m}{dk} \right \rfloor} j]]$

可以发现$\sum_{i=1}^{\left \lfloor \frac{n}{dk} \right \rfloor} i$和$\sum_{j=1}^{\left \lfloor \frac{m}{dk} \right \rfloor} j$都是等差数列求和，$\sum_{k=1}^{\left \lfloor \frac{n}{d} \right \rfloor} k^{2}\mu(k)$可以通过预处理前缀和来求，两次整数分块即可

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
57
58
59
60
61
62
63
64
65
66
67
#include <iostream>
#include <algorithm>
#include <cstring>
#include <cstdio>

using namespace std;

const int N = 10000010;
const int mod = 20101009;

int isp[N], p[N], sum[N], mu[N];
int tot, n, m;

void init(int n)
{
    mu[1] = 1;
    for (int i = 2; i <= n; i++) {
        if (!isp[i]) {
            p[++tot] = i;
            mu[i] = -1;
        }
        for (int j = 1; j <= tot && p[j] <= n / i; j++) {
            isp[i * p[j]] = 1;
            if (0 == i % p[j]) {
                mu[i * p[j]] = 0;
                break;
            }
            mu[i * p[j]] = -mu[i];
        }
    }
    for (int i = 1; i <= n; i++)
        sum[i] = (sum[i - 1] + 1ll * i * i % mod * (mu[i] + mod)) % mod;
}

int qs(int x, int y)
{
    return (1ll * x * (x + 1) / 2 % mod) * (1ll * y * (y + 1) / 2 % mod) % mod;
}

int calc(int n, int m)
{
    int res = 0;
    for (int l = 1, r = 1; l <= n; l = r + 1) {
        r = min(n / (n / l), m / (m / l));
        res = (res + 1ll * (sum[r] - sum[l - 1] + mod) * qs(n / l, m / l) % mod) % mod;
    }
    return res;
}

int solve(int n, int m)
{
    int res = 0;
    for (int l = 1, r = 1; l <= n; l = r + 1) {
        r = min(n / (n / l), m / (m / l));
        res = (res + 1ll * (r - l + 1) * (l + r) / 2 % mod * calc(n / l, m / l) % mod) % mod;
    }
    return res;
}

int main()
{
    init(N - 5);
    scanf("%d%d", &n, &m);
    if (n > m) swap(n, m);
    printf("%d\n", solve(n, m));
    return 0;
}

AT5200 [AGC038C] LCMs

题目链接：AtCoder - agc038_c

题意：给定一个长度为n得序列$A_{1}，A_{2}，…，A_{n}$，求$\sum_{i=1}^{n} \sum_{j=i+1}^{n} lcm(A_{i},A_{j})$

思路：令M=max($A_{i}$)，$\sum_{i=1}^{n} \sum_{j=i+1}^{n} lcm(A_{i},A_{j})=\frac{\sum_{i=1}^{n} \sum_{j=1}^{n} lcm(A_{i},A_{j})-\sum_{i=1}^{n}A_{i}}{2}$

考虑如何求$res=\sum_{i=1}^{n} \sum_{j=1}^{n} lcm(A_{i},A_{j})$

$\sum_{i=1}^{n} \sum_{j=1}^{n} lcm(A_{i},A_{j})=\sum_{i=1}^{n} \sum_{j=1}^{n} \frac{A_{i}A_{j}}{gcd(A_{i},A_{j})}$

令$d=gcd(A_{i},A_{j})$得$\sum_{i=1}^{n} \sum_{j=1}^{n} lcm(A_{i},A_{j})=\sum_{d=1}^{M} \frac{1}{d} \sum_{i=1}^{n} \sum_{j=1}^{n} A_{i}A_{j}[gcd(A_{i},A_{j})=d]=$

$\sum_{d=1}^{M} \frac{1}{d} \sum_{i=1}^{n} \sum_{j=1}^{n} A_{i}A_{j}[gcd(\frac{A_{i}}{d},\frac{A_{j}}{d})=1][d|A_{i}][d|A_{j}]$

根据莫比乌斯反演的性质得$res=\sum_{d=1}^{M} \frac{1}{d} \sum_{i=1}^{n} \sum_{j=1}^{n} A_{i}A_{j} \sum_{k|gcd(\frac{A_{i}}{d},\frac{A_{j}}{d})} \mu(k) [d|A_{i}][d|A_{j}]$

将k提到前面得$res=\sum_{d=1}^{M} \frac{1}{d} \sum_{k=1}^{\left \lfloor \frac{M}{d} \right \rfloor} \mu(k) \sum_{i=1}^{n} \sum_{j=1}^{n} A_{i}A_{j} [kd|A_{i}][kd|A_{j}]$

所以最后$res=\sum_{d=1}^{M} \frac{1}{d} \sum_{k=1}^{\left \lfloor \frac{M}{d} \right \rfloor} \mu(k) (\sum_{i=1}^{n}[kd|A_{i}]A_{i})^{2}$

设$f(t)=(\sum_{i=1}^{n}[t|A_{i}]A_{i})^{2}$，显然$f(t)$可以在$O(nlnn)$的时间内预处理出来，再通过枚举d和k暴力求解即可

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
57
58
59
60
61
#include <iostream>
#include <algorithm>
#include <cstring>
#include <cstdio>

using namespace std;

typedef long long ll;

const int N = 1000010;
const ll mod = 998244353;

int isp[N], p[N], mu[N], n, tot, a[N], imax;
ll inv[N], f[N], c[N], sum;

void init(int n)
{
    mu[1] = 1;
    for (int i = 2; i <= n; i++) {
        if (!isp[i]) {
            p[++tot] = i;
            mu[i] = -1;
        }
        for (int j = 1; j <= tot && p[j] <= n / i; j++) {
            isp[i * p[j]] = 1;
            if (0 == i % p[j]) {
                mu[i * p[j]] = 0;
                break;
            }
            else mu[i * p[j]] = -mu[i];
        }
    }
    inv[0] = inv[1] = 1;
    for (int i = 2; i <= n; i++) inv[i] = (mod - mod / i) * inv[mod % i] % mod;
}

int main()
{
    scanf("%d", &n);
    for (int i = 1; i <= n; i++) {
        scanf("%lld", &a[i]);
        imax = max(imax, a[i]);
        sum = (sum + a[i]) % mod;
        c[a[i]] += 1;
    }
    init(imax);
    for (int i = 1; i <= imax; i++) {
        for (int j = i; j <= imax; j += i)
            f[i] = (f[i] + c[j] * j) % mod;
        f[i] = (f[i] * f[i]) % mod;
    }
    ll res = 0;
    for (int d = 1; d <= imax; d++) {
        for (int k = 1, t = d; t <= imax; t += d, k++) {
            res = (res + inv[d] * mu[k] % mod * f[t] % mod) % mod;
        }
    }
    res = ((res - sum) % mod + mod) % mod * inv[2] % mod;
    printf("%lld\n", res);
    return 0;
}