「学习笔记」莫比乌斯反演

发表于 更新于 分类于 阅读次数： 本文字数： 14k 阅读时长 ≈ 28 分钟

μ∗1=e\mu * 1=e

φ∗1=I\varphi * 1=I

Luogu P4450 双亲数

Description

给定 n,m,dn,m,d，qiu

∑i=1n∑j=1m[gcd⁡(i,j)=d]\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=d]

1≤n,m≤106,1≤d≤min(n,m)1\le n,m\le 10^6,1\le d\le min(n,m)

Solution

∑i=1n∑j=1m[gcd⁡(i,j)=d]=∑i=1nd∑j=1md[gcd⁡(i,j)=1]=∑i=1nd∑j=1md∑k∣gcd(i,j)μ(k)=∑k=1n∑i=1nkd∑j=1mkdμ(k)=∑k=1n⌊nkd⌋⌊mkd⌋μ(k)\begin{aligned} \sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=d] &=\sum_{i=1}^{\frac n d}\sum_{j=1}^{\frac m d}[\gcd(i,j)=1] \\ &=\sum_{i=1}^{\frac n d}\sum_{j=1}^{\frac m d}\sum_{k|gcd(i,j)}\mu(k) \\ &=\sum_{k=1}^n\sum_{i=1}^{\frac n {kd}}\sum_{j=1}^{\frac m {kd}}\mu(k) \\ &=\sum_{k=1}^n\left\lfloor \frac n {kd} \right\rfloor \left\lfloor \frac m {kd} \right\rfloor \mu(k) \end{aligned}

然后就可以整除分块了

Code 
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
#include <bits/stdc++.h>
#define ll long long
#define db double
#define gc getchar
#define pc putchar

using namespace std;

namespace IO
{
    template <typename T>
    void read(T &x)
    {
        x = 0; bool f = 0; char c = gc();
        while(!isdigit(c)) f |= c == '-', c = gc();
        while(isdigit(c)) x = x * 10 + c - '0', c = gc();
        if(f) x = -x;
    }

    template <typename T>
    void write(T x)
    {
        if(x < 0) pc('-'), x = -x;
        if(x > 9) write(x / 10);
        pc('0' + x % 10);
    }
}
using namespace IO;

const int N = 1e6 + 5;

int p[N], tot, mu[N];
bool flag[N];

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

ll solve(int n, int m, int d)
{
    n /= d, m /= d;
    if(n > m) n ^= m ^= n ^= m;
    ll res = 0;
    for(int l = 1, r; l <= n; l = r + 1)
    {
        r = min(n / (n / l), m / (m / l));
        res += 1ll * (mu[r] - mu[l - 1]) * (n / l) * (m / l);
    }
    return res;
}

int main()
{
    init(1e6);
    int n, m, d;;
    read(n), read(m), read(d);
    write(solve(n, m, d)), pc('\n');
    return 0;
}
// A.S.

还有一道双倍经验 Luogu P3455 [POI2007]ZAP-Queries

Luogu P2522 Problem b

Description

nn 组询问，每次询问给定 a,b,c,d,ka,b,c,d,k，求

∑i=ab∑j=cd[gcd⁡(i,j)=k]\sum_{i=a}^b\sum_{j=c}^d[\gcd(i,j)=k]

1≤n,k≤5×104,1≤a≤b≤5×104,1≤c≤d≤5×1041\le n,k \le 5\times 10^4,1\le a\le b\le 5\times 10^4,1\le c\le d\le 5\times 10^4

Solution

不难想到容斥，于是题目转化为求

∑i=1n∑j=1m[gcd⁡(i,j)=k]\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=k]

就跟上一题一样了

Code 
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
#include <bits/stdc++.h>
#define ll long long
#define db double
#define gc getchar
#define pc putchar

using namespace std;

namespace IO
{
    template <typename T>
    void read(T &x)
    {
        x = 0; bool f = 0; char c = gc();
        while(!isdigit(c)) f |= c == '-', c = gc();
        while(isdigit(c)) x = x * 10 + c - '0', c = gc();
        if(f) x = -x;
    }

    template <typename T>
    void write(T x)
    {
        if(x < 0) pc('-'), x = -x;
        if(x > 9) write(x / 10);
        pc('0' + x % 10);
    }
}
using namespace IO;

const int N = 5e4 + 5;
int p[N], tot, mu[N], sum[N];
bool flag[N];

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

int k;

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

int main()
{
    init(N - 5);
    int T; read(T);
    while(T--)
    {
        int a, b, c, d;
        read(a), read(b), read(c), read(d), read(k);
        write(calc(b, d) - calc(a - 1, d) - calc(b, c - 1) + calc(a - 1, c - 1)), pc('\n');
    }
    return 0;
}
// A.S.

Luogu P3911 最小公倍数之和

Description

给定 a1,a2,…ana_1,a_2,\dots a_n，求

∑i=1n∑j=1nlcm(ai,aj)\sum_{i=1}^n\sum_{j=1}^n\text{lcm}(a_i,a_j)

1≤n≤5×104,1≤ai≤5×1041\le n \le 5\times 10^4,1\le a_i \le 5 \times 10^4

Solution

设 m=max(ai)m=max(a_i)，cic_i 表示 ∑a\sum a 中有多少个 ii。

∑i=1n∑j=1nlcm(ai,aj)=∑i=1n∑j=1nai⋅ajgcd⁡(ai,aj)=∑i=1m∑j=1mi⋅jgcd⁡(i,j)⋅ci⋅cj=∑d=1m∑i=1m∑j=1m[gcd⁡(i,j)=d]i⋅jd⋅ci⋅cj=∑d=1m∑i=1md∑j=1md[gcd⁡(i,j)=1]i⋅j⋅d⋅cid⋅cjd=∑d=1md∑i=1md∑j=1md∑k∣gcd⁡(i,j)μ(k)⋅i⋅j⋅cidk⋅cjdk=∑d=1md∑k=1mdμ(k)⋅k2∑i=1mdk∑j=1mdki⋅j⋅cidk⋅cjdk=∑dk=1mdk∑k∣dkμ(k)⋅k∑i=1mdk(i⋅cidk)2\begin{aligned} \sum_{i=1}^n\sum_{j=1}^n\text{lcm}(a_i,a_j) &=\sum_{i=1}^n\sum_{j=1}^n\dfrac{a_i·a_j}{\gcd(a_i,a_j)} \\ &=\sum_{i=1}^m\sum_{j=1}^m\dfrac{i·j}{\gcd(i,j)}·c_i·c_j \\ &=\sum_{d=1}^m\sum_{i=1}^m\sum_{j=1}^m[\gcd(i,j)=d]\dfrac{i·j}{d}·c_i·c_j \\ &=\sum_{d=1}^m\sum_{i=1}^{\frac m d}\sum_{j=1}^{\frac m d}[\gcd(i,j)=1]i·j·d·c_{id}·c_{jd} \\ &=\sum_{d=1}^md\sum_{i=1}^{\frac m d}\sum_{j=1}^{\frac m d}\sum_{k|\gcd(i,j)}\mu(k)·i·j·c_{idk}·c_{jdk} \\ &=\sum_{d=1}^md\sum_{k=1}^{\frac m d}\mu(k)·k^2\sum_{i=1}^{\frac m{dk}}\sum_{j=1}^{\frac m{dk}}i·j·c_{idk}·c_{jdk} \\ &=\sum_{dk=1}^mdk\sum_{k|dk}\mu(k)·k\sum_{i=1}^{\frac m{dk}}(i·c_{idk})^2 \\ \end{aligned}

令 T=dkT=dk

∑i=1n∑j=1nlcm(ai,aj)=∑T=1mT∑k∣Tμ(k)⋅k∑i=1mT(i⋅ciT)2\sum_{i=1}^n\sum_{j=1}^n\text{lcm}(a_i,a_j)= \sum_{T=1}^mT\sum_{k|T}\mu(k)·k\sum_{i=1}^{\frac m T}(i·c_{iT})^2 \\

然后第二、三个求和可以预处理，复杂度是调和级数，最后再 O(m)O(m) 统计答案即可。

Code 
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
#include <bits/stdc++.h>
#define ll long long
#define db double
#define gc getchar
#define pc putchar

using namespace std;

namespace IO
{
    template <typename T>
    void read(T &x)
    {
        x = 0; bool f = 0; char c = gc();
        while(!isdigit(c)) f |= c == '-', c = gc();
        while(isdigit(c)) x = x * 10 + c - '0', c = gc();
        if(f) x = -x;
    }

    template <typename T>
    void write(T x)
    {
        if(x < 0) pc('-'), x = -x;
        if(x > 9) write(x / 10);
        pc('0' + x % 10);
    }
}
using namespace IO;

const int N = 5e4 + 5;

int n, m, a[N];
int p[N], tot, mu[N];
bool flag[N];
ll c[N], sum[N], smu[N];

void init()
{
    mu[1] = 1;
    for(int i = 2; i <= m; i++)
    {
        if(!flag[i]) p[++tot] = i, mu[i] = -1;
        for(int j = 1; j <= tot && i * p[j] <= m; j++)
        {
            flag[i * p[j]] = 1;
            if(i % p[j] == 0)
            {
                mu[i * p[j]] = 0;
                break;
            }
            mu[i * p[j]] = -mu[i];
        }
    }
    for(int i = 1; i <= m; i++)
        for(int j = i; j <= m; j += i)
            smu[j] = smu[j] + mu[i] * i;
    for(int i = 1; i <= n; i++) c[a[i]]++; 
    for(int i = 1; i <= m; i++)
        for(int j = 1; j <= m / i; j++)
            sum[i] = sum[i] + j * c[i * j];
    for(int i = 1; i <= m; i++) sum[i] = sum[i] * sum[i];
    return;
}

signed main()
{
    read(n);
    for(int i = 1; i <= n; i++) read(a[i]), m = max(m, a[i]);
    init();
    ll ans = 0;
    for(int i = 1; i <= m; i++)
        ans += i * smu[i] * sum[i];
    write(ans), pc('\n');
    return 0;
}
// A.S.

AT5200 LCMs

Description

给定 a1,a2,…ana_1,a_2,\dots a_n，求

∑i=1n∑j=i+1nlcm(ai,aj)\sum_{i=1}^n\sum_{j=i+1}^n\text{lcm}(a_i,a_j)

对 998244353998244353 取模

1≤n≤2×105,1≤ai≤1051\le n \le 2\times 10^5,1\le a_i \le 10^5

Solution

∑i=1n∑j=i+1nlcm(ai,aj)=∑i=1n∑j=1nlcm(ai,aj)−sumi=1nai2\sum_{i=1}^n\sum_{j=i+1}^n\text{lcm}(a_i,a_j) =\dfrac{\sum_{i=1}^n\sum_{j=1}^n\text{lcm}(a_i,a_j)-sum_{i=1}^na_i}{2}

剩下的就跟上一题一样了。

注意有 μ\mu 可能会有负数，所以还是老老实实取模吧，不然看不出来要调好久 qaq

COde 
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
#include <bits/stdc++.h>
#define ll long long
#define db double
#define gc getchar
#define pc putchar

using namespace std;

namespace IO
{
    template <typename T>
    void read(T &x)
    {
        x = 0; bool f = 0; char c = gc();
        while(!isdigit(c)) f |= c == '-', c = gc();
        while(isdigit(c)) x = x * 10 + c - '0', c = gc();
        if(f) x = -x;
    }

    template <typename T>
    void write(T x)
    {
        if(x < 0) pc('-'), x = -x;
        if(x > 9) write(x / 10);
        pc('0' + x % 10);
    }
}
using namespace IO;

const int N = 2e5 + 5;
const int M = 1e6 + 5;
const int mod = 998244353;

ll add(ll x) {return x < mod ? x : x - mod;}
ll sub(ll x) {return x < 0 ? x : x + mod;}
ll qpow(ll a, int b)
{
    ll res = 1;
    while(b)
    {
        if(b & 1) res = res * a % mod;
        a = a * a % mod, b >>= 1;
    }
    return res;
}

int n, m, a[N];
int p[M], tot, mu[M];
bool flag[M];
ll c[M], sum[M], smu[M];

void init()
{
    mu[1] = 1;
    for(int i = 2; i <= m; i++)
    {
        if(!flag[i]) p[++tot] = i, mu[i] = -1;
        for(int j = 1; j <= tot && i * p[j] <= m; j++)
        {
            flag[i * p[j]] = 1;
            if(i % p[j] == 0)
            {
                mu[i * p[j]] = 0;
                break;
            }
            mu[i * p[j]] = -mu[i];
        }
    }
    for(int i = 1; i <= m; i++)
        for(int j = i; j <= m; j += i)
            smu[j] = add((smu[j] + mu[i] * i % mod) % mod + mod);
    for(int i = 1; i <= n; i++) c[a[i]]++; 
    for(int i = 1; i <= m; i++)
        for(int j = 1; j <= m / i; j++)
            sum[i] = add(sum[i] + j * c[i * j] % mod);
    for(int i = 1; i <= m; i++) sum[i] = sum[i] * sum[i] % mod;
    return;
}

signed main()
{
    read(n);
    for(int i = 1; i <= n; i++) read(a[i]), m = max(m, a[i]);
    init();
    ll ans = 0;
    for(int i = 1; i <= m; i++)
        ans = ((ans + i * smu[i] % mod * sum[i] % mod) % mod + mod) % mod;
    for(int i = 1; i <= n; i++) ans = ((ans - a[i]) % mod + mod) % mod;
    ans = ans * qpow(2, mod - 2) % mod;
    write(ans), pc('\n');
    return 0;
}
// A.S.

Luogu P5221 Product

Description

给定 nn，求

∏i=1n∏j=1nlcm(i,j)gcd⁡(i,j)( mod  104857601)\prod_{i=1}^n\prod_{j=1}^n\dfrac{\text{lcm}(i,j)}{\gcd(i,j)}\quad(\bmod\ 104857601)

1≤n≤1061\le n \le 10^6

Solution

∏i=1n∏j=1nlcm(i,j)gcd⁡(i,j)=∏i=1n∏j=1ni∗jgcd⁡(i,j)2=(n!)2∏i=1n∏j=1ngcd⁡(i,j)2∏i=1n∏j=1ngcd⁡(i,j)2=(∏i=1n∏j=1ngcd⁡(i,j))2∏i=1n∏j=1ngcd⁡(i,j)=∏d=1n∏i=1n∏j=1n[gcd⁡(i,j)=d]d=∏d=1nd∑i=1n∑j=1n[gcd⁡(i,j)=d]∑i=1n∑j=1n[gcd⁡(i,j)=d]=∑i=1nd∑j=1nd[gcd⁡(i,j)=1]=1+2∑i=1nφ(i)\begin{aligned} \prod_{i=1}^n\prod_{j=1}^n\dfrac{\text{lcm}(i,j)}{\gcd(i,j)} &=\prod_{i=1}^n\prod_{j=1}^n\dfrac{i*j}{\gcd(i,j)^2} \\ &=\dfrac{(n!)^2}{\prod_{i=1}^n\prod_{j=1}^n\gcd(i,j)^2} \\ \prod_{i=1}^n\prod_{j=1}^n\gcd(i,j)^2 &=(\prod_{i=1}^n\prod_{j=1}^n\gcd(i,j))^2 \\ \prod_{i=1}^n\prod_{j=1}^n\gcd(i,j) &=\prod_{d=1}^n\prod_{i=1}^n\prod_{j=1}^n[\gcd(i,j)=d]d \\ &=\prod_{d=1}^nd^{\sum_{i=1}^n\sum_{j=1}^n[\gcd(i,j)=d]} \\ \sum_{i=1}^n\sum_{j=1}^n[\gcd(i,j)=d] &=\sum_{i=1}^{\frac n d}\sum_{j=1}^{\frac n d}[\gcd(i,j)=1] \\ &=1+2\sum_{i=1}^n\varphi(i) \end{aligned}

预处理一下 φ\varphi 的前缀和即可，设 sum(i)=∑j=1iφ(j)sum(i)=\sum_{j=1}^i\varphi(j)。

∏i=1n∏j=1nlcm(i,j)gcd⁡(i,j)=(n!)2(∏d=1nd2∗sum(nd)+1)2\prod_{i=1}^n\prod_{j=1}^n\dfrac{\text{lcm}(i,j)}{\gcd(i,j)} =\dfrac{(n!)^2}{(\prod_{d=1}^n d^{2*sum(\frac n d)+1})^2}

根据欧拉定理，可以将 2∗sum(nd)+1 mod  (p−1)2*sum(\frac n d)+1\bmod\ (p-1)

Code 
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
#include <bits/stdc++.h>
#define ll long long
#define db double
#define gc getchar
#define pc putchar

using namespace std;

namespace IO
{
    template <typename T>
    void read(T &x)
    {
        x = 0; bool f = 0; char c = gc();
        while(!isdigit(c)) f |= c == '-', c = gc();
        while(isdigit(c)) x = x * 10 + c - '0', c = gc();
        if(f) x = -x;
    }

    template <typename T>
    void write(T x)
    {
        if(x < 0) pc('-'), x = -x;
        if(x > 9) write(x / 10);
        pc('0' + x % 10);
    }
}
using namespace IO;

const int N = 1e6 + 5;
const int p = 104857601;

ll qpow(ll a, int b)
{
    ll res = 1;
    while(b)
    {
        if(b & 1) res = res * a % p;
        a = a * a % p, b >>= 1;
    }
    return res;
}

int pri[N], tot, phi[N];
bool flag[N];

void getphi(int n)
{
    phi[1] = 1;
    for(int i = 2; i <= n; i++)
    {
        if(!flag[i]) pri[++tot] = i, phi[i] = i - 1;
        for(int j = 1; j <= tot && i * pri[j] <= n; j++)
        {
            flag[i * pri[j]] = 1;
            if(i % pri[j] == 0)
            {
                phi[i * pri[j]] = phi[i] * pri[j];
                break;
            }
            phi[i * pri[j]] = phi[i] * (pri[j] - 1);
        }
    }
    for(int i = 1; i <= n; i++) phi[i] = (phi[i - 1] + phi[i]) % (p - 1);
    return;
}

int main()
{
    int n; read(n);
    getphi(n);
    int ans = 1;
    for(int i = 1; i <= n; i++) ans = 1ll * ans * i % p;
    ans = qpow(ans, 2 * n);
    int tmp = 1;
    for(int i = 1; i <= n; i++)
        tmp = 1ll * tmp * qpow(i, (2 * phi[n / i] - 1) % (p - 1)) % p;
    ans = 1ll * ans * qpow(1ll * tmp * tmp % p, p - 2) % p;
    write(ans), pc('\n');
    return 0;
}
// A.S.

Luogu P4917 天守阁的地板

Description

TT 组数据，每组数据给定一个 nn，求

∏i=1n∏j=1nlcm(i,j)gcd⁡(i,j)( mod  19260817)\prod_{i=1}^n\prod_{j=1}^n\dfrac{\text{lcm}(i,j)}{\gcd(i,j)}\quad(\bmod\ 19260817)

1≤T≤1000,1≤n≤1061\le T \le 1000,1\le n \le 10^6

有没有发现跟上一题一样（强烈谴责题面

Solution

本题与上一题唯一的区别就是最后计算答案时需要整除分块，因为 O(Tn)O(Tn) 过不去，而我上一题偷懒没写（

Code 
1
2
3
4
5
6
7
8
9
10
int calc(int n)
{
    int ans = qpow(fac[n], n << 1), res = 1;
    for(int l = 1, r; l <= n; l = r + 1)
    {
        r = n / (n / l);
        res = 1ll * res * qpow(1ll * fac[r] * qpow(fac[l - 1], p - 2) % p, 2 * phi[n / l] - 1) % p;
    }
    return 1ll * ans * qpow(1ll * res * res % p, p - 2) % p;
}

Luogu P3327 约数个数和

Description

设 d(x)d(x) 为 xx 的约数个数，给定 n,mn,m，求

∑i=1n∑j=1md(i∗j)\sum_{i=1}^n\sum_{j=1}^md(i*j)

TT 组数据

1≤T,n,m≤500001\le T,n,m\le 50000

Solution

首先有个性质

d(x∗y)=∑i∣x∑j∣y[gcd⁡(i,j)=1]d(x*y)=\sum_{i|x}\sum_{j|y}[\gcd(i,j)=1]

所以原式为

∑x=1n∑y=1m∑i∣x∑j∣y[gcd⁡(i,j)=1]=∑i=1n∑j=1m⌊ni⌋⌊mj⌋[gcd⁡(i,j)=1]\sum_{x=1}^n\sum_{y=1}^m\sum_{i|x}\sum_{j|y}[\gcd(i,j)=1] =\sum_{i=1}^n\sum_{j=1}^m\left\lfloor\dfrac n i\right\rfloor \left\lfloor\dfrac m j\right\rfloor[\gcd(i,j)=1]

设

f(x)=∑i=1n∑j=1m⌊ni⌋⌊mj⌋[gcd⁡(i,j)=x]g(x)=∑x∣df(d)f(x)=\sum_{i=1}^n\sum_{j=1}^m\left\lfloor\dfrac n i\right\rfloor\left\lfloor\dfrac m j\right\rfloor[\gcd(i,j)=x] \\ g(x)=\sum_{x|d}f(d)

那么

g(x)=∑i=1n∑j=1m⌊ni⌋⌊mj⌋[x∣gcd⁡(i,j)]=∑i=1nx∑j=1mx⌊nix⌋⌊mjx⌋\begin{aligned} g(x)&=\sum_{i=1}^n\sum_{j=1}^m\left\lfloor\dfrac n i\right\rfloor\left\lfloor\dfrac m j\right\rfloor[x|\gcd(i,j)] \\ &=\sum_{i=1}^{\frac n x}\sum_{j=1}^{\frac m x}\left\lfloor\dfrac n {ix}\right\rfloor\left\lfloor\dfrac m {jx}\right\rfloor \end{aligned}

因为

f(n)=∑n∣dμ(nd)g(d)f(n)=\sum_{n|d}\mu(\frac n d)g(d)

所以

f(1)=∑1∣dμ(d)g(d)=∑i=1nμ(i)g(i)\begin{aligned} f(1)&=\sum_{1|d}\mu(d)g(d) \\ &=\sum_{i=1}^n\mu(i)g(i) \end{aligned}

考虑如何求 g(x)g(x)，设 s(x)=∑i=1n⌊ni⌋s(x)=\sum_{i=1}^n\lfloor\frac n i\rfloor

g(x)=s(⌊ni⌋)∗s(⌊mj⌋)g(x)=s(\lfloor\frac n i\rfloor)*s(\lfloor\frac m j\rfloor)

Code 
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
#include <bits/stdc++.h>
#define ll long long
#define db double
#define gc getchar
#define pc putchar

using namespace std;

namespace IO
{
    template <typename T>
    void read(T &x)
    {
        x = 0; bool f = 0; char c = gc();
        while(!isdigit(c)) f |= c == '-', c = gc();
        while(isdigit(c)) x = x * 10 + c - '0', c = gc();
        if(f) x = -x;
    }

    template <typename T>
    void write(T x)
    {
        if(x < 0) pc('-'), x = -x;
        if(x > 9) write(x / 10);
        pc('0' + x % 10);
    }
}
using namespace IO;

const int N = 5e4 + 5;

int p[N], tot;
ll mu[N], sum[N];
bool flag[N];

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

    for(int i = 1; i <= n; i++)
    {
        for(int l = 1, r; l <= i; l = r + 1)
        {
            r = i / (i / l);
            sum[i] += 1ll * (r - l + 1) * (i / l);
        }
    }
    return;
}

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

int main()
{
    init(5e4);
    int T; read(T);
    while(T--)
    {
        int n, m;
        read(n), read(m);
        write(calc(n, m)), pc('\n');
    }
    return 0;
}
// A.S.

Luogu P1829 Crash的数字表格 / JZPTAB

Description

给定 n,mn,m，求

∑i=1n∑j=1mlcm(i,j)( mod  20101009)\sum_{i=1}^n\sum_{j=1}^m\text{lcm}(i,j)\quad(\bmod\ 20101009)

1≤n,m≤1071\le n,m\le 10^7

Solution

默认 n≤mn\le m，除法向下取整

∑i=1n∑j=1mlcm(i,j)=∑i=1n∑j=1mi∗jgcd⁡(i,j)=∑i=1n∑j=1m∑d=1n[gcd⁡(i,j)=d]i∗jd=∑d=1nd∑i=1nd∑j=1md[gcd⁡(i,j)=1]i∗j\begin{aligned} \sum_{i=1}^n\sum_{j=1}^m\text{lcm}(i,j) &=\sum_{i=1}^n\sum_{j=1}^m\dfrac{i*j}{\gcd(i,j)} \\ &=\sum_{i=1}^n\sum_{j=1}^m\sum_{d=1}^{n}[\gcd(i,j)=d]\dfrac{i*j}{d}\\ &=\sum_{d=1}^nd\sum_{i=1}^{\frac n d}\sum_{j=1}^{\frac m d}[\gcd(i,j)=1]i*j \\ \end{aligned}

设

sum(n,m)=∑i=1n∑j=1m[gcd⁡(i,j)=1]i∗j=∑i=1n∑j=1m∑d∣gcd⁡(i,j)μ(d)∗i∗j=∑d=1nμ(d)∗d2∑i=1nd∑j=1mdi∗j\begin{aligned} sum(n,m)&=\sum_{i=1}^{n}\sum_{j=1}^{m}[\gcd(i,j)=1]i*j \\ &=\sum_{i=1}^n\sum_{j=1}^m\sum_{d|\gcd(i,j)}\mu(d)*i*j \\ &=\sum_{d=1}^n\mu(d)*d^2\sum_{i=1}^{\frac n d}\sum_{j=1}^{\frac m d}i*j \\ \end{aligned}

前面预处理前缀和即可，后面就是两个等差数列的和的乘积

记

getsum(n,m)=n(n+1)2×m(m+1)2getsum(n,m)=\dfrac{n(n+1)} 2\times \dfrac{m(m+1)}2

则

sum(n,m)=∑d=1nμ(d)∗d2∗getsum(nd,md)sum(n,m)=\sum_{d=1}^n\mu(d)*d^2*getsum(\frac n d,\frac m d)

可以整除分块，然后原式为

∑d=1nd∗sum(nd,md)\sum_{d=1}^nd*sum(\frac n d,\frac m d)

发现又可以整除分块

Code 
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
#include <bits/stdc++.h>
#define int long long
#define db double
#define gc getchar
#define pc putchar

using namespace std;

namespace IO
{
    template <typename T>
    void read(T &x)
    {
        x = 0; bool f = 0; char c = gc();
        while(!isdigit(c)) f |= c == '-', c = gc();
        while(isdigit(c)) x = x * 10 + c - '0', c = gc();
        if(f) x = -x;
    }

    template <typename T>
    void write(T x)
    {
        if(x < 0) pc('-'), x = -x;
        if(x > 9) write(x / 10);
        pc('0' + x % 10);
    }
}
using namespace IO;

const int N = 1e7 + 5;
const int mod = 20101009;

int p[N], tot, mu[N];
int sum[N];
bool flag[N];

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

int getsum(int n, int m)
{
    return (n * (n + 1) / 2 % mod) * (m * (m + 1) / 2 % mod) % mod;
}

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

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

signed main()
{
    int n, m;
    read(n), read(m);
    if(n > m) n ^= m ^= n ^= m;
    init(n, m);
    write(solve(n, m)), pc('\n');
    return 0;
}
// A.S.

Luogu P3768 简单的数学题

Description

给定 n,pn,p，求

(∑i=1n∑j=1nijgcd⁡(i,j)) mod  p\left( \sum_{i=1}^n\sum_{j=1}^nij\gcd(i,j) \right)\bmod\ p

n≤1010,5×108≤p≤1.1×109n\le 10^{10},5\times 10^8\le p \le 1.1\times 10^9

Solution

∑i=1n∑j=1nijgcd⁡(i,j)=∑d=1nd∑i=1n∑j=1n[gcd⁡(i,j)=d]ij=∑d=1nd3∑i=1nd∑j=1nd[gcd⁡(i,j)=1]ij=∑d=1nd3∑k=1ndμ(k)k2(∑i=1ndki)2\begin{aligned} \sum_{i=1}^n\sum_{j=1}^nij\gcd(i,j) &=\sum_{d=1}^nd\sum_{i=1}^n\sum_{j=1}^n[\gcd(i,j)=d]ij \\ &=\sum_{d=1}^nd^3\sum_{i=1}^{\frac n d}\sum_{j=1}^{\frac n d}[\gcd(i,j)=1]ij \\ &=\sum_{d=1}^nd^3\sum_{k=1}^{\frac n d}\mu(k)k^2(\sum_{i=1}^{\frac n {dk}}i)^2 \end{aligned}

设 sum(n)=∑i=1nisum(n)=\sum_{i=1}^ni

令 T=dkT=dk

∑T=1nsum(nT)2∑d∣Td3(Td)2μ(Td)∑T=1nsum(nT)2T2∑d∣Tdμ(Td)\sum_{T=1}^nsum(\frac n T)^2\sum_{d|T}d^3(\frac T d)^2\mu(\frac T d) \\ \sum_{T=1}^nsum(\frac n T)^2T^2\sum_{d|T}d\mu(\frac T d)

现在只需要筛出来后面的那一部分

id∗μ=φ∑d∣Tdμ(Td)=φ(T)T2∑d∣Tdμ(Td)=T2φ(T)id*\mu=\varphi \\ \sum_{d|T}d\mu(\frac T d)=\varphi(T) \\ T^2\sum_{d|T}d\mu(\frac T d)=T^2\varphi(T)

设

f(i)=i2φ(i),S(n)=∑i=1nf(i)f(i)=i^2\varphi(i),S(n)=\sum_{i=1}^nf(i)

看一下杜教筛套路式子

g(1)S(n)=∑i=1n(g∗f)(i)−∑i=2ng(i)S(ni)g(1)S(n)=\sum_{i=1}^n(g*f)(i)-\sum_{i=2}^ng(i)S(\frac n i)

考虑 gg 应该是什么

因为 φ∗1=id\varphi*1=id，令 g(i)=i2g(i)=i^2

S(n)=∑i=1n(g∗f)(i)−∑i=2ng(i)S(ni)(g∗f)(i)=∑d∣if(d)g(id)=∑d∣id2φ(d)(id)2=i2∑d∣iφ(d)=i3S(n)=\sum_{i=1}^n(g*f)(i)-\sum_{i=2}^ng(i)S(\frac n i) \\ \begin{aligned} (g*f)(i)&=\sum_{d|i}f(d)g(\frac i d) \\ &=\sum_{d|i}d^2\varphi(d)(\frac i d)^2 \\ &=i^2\sum_{d|i}\varphi(d) \\ &=i^3 \end{aligned}

所以

S(n)=∑i=1ni3−∑i=2ni2S(ni)S(n)=\sum_{i=1}^ni^3-\sum_{i=2}^ni^2S(\frac n i)

其中

∑i=1ni3=(∑i=1ni)2\sum_{i=1}^ni^3=(\sum_{i=1}^ni)^2

因此原式

∑T=1nsum(nT)2f(T)\sum_{T=1}^nsum(\frac n T)^2f(T)

就可以整除分块加杜教筛解决了

Code 
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
105
106
107
108
#include <bits/stdc++.h>
#define ll long long
#define db double
#define gc getchar
#define pc putchar

using namespace std;

namespace IO
{
    template <typename T>
    void read(T &x)
    {
        x = 0; bool f = 0; char c = gc();
        while(!isdigit(c)) f |= c == '-', c = gc();
        while(isdigit(c)) x = x * 10 + c - '0', c = gc();
        if(f) x = -x;
    }

    template <typename T>
    void write(T x)
    {
        if(x < 0) pc('-'), x = -x;
        if(x > 9) write(x / 10);
        pc('0' + x % 10);
    }
}
using namespace IO;

const int N = 4e6 + 5;

ll p, n;
ll pri[N], tot, phi[N];
ll sum[N], inv6;
bool flag[N];
map <ll, ll> mp;

ll add(ll x) {return x < p ? x : x - p;}
ll sub(ll x) {return x < 0 ? x + p : x;}
ll s2(ll x) {x %= p; return x * (x + 1) % p * (2 * x + 1) % p * inv6 % p;}
ll s3(ll x) {x %= p; return (x * (x + 1) / 2 % p) * (x * (x + 1) / 2 % p) % p;}

ll qpow(ll a, int b)
{
    ll res = 1;
    while(b)
    {
        if(b & 1) res = res * a % p;
        a = a * a % p, b >>= 1;
    }
    return res;
}

void init(int n)
{
    phi[1] = 1;
    for(int i = 2; i < n; i++)
    {
        if(!flag[i]) pri[++tot] = i, phi[i] = i - 1;
        for(int j = 1; j <= tot && i * pri[j] < n; j++)
        {
            flag[i * pri[j]] = 1;
            if(i % pri[j] == 0)
            {
                phi[i * pri[j]] = phi[i] * pri[j];
                break;
            }
            phi[i * pri[j]] = phi[i] * (pri[j] - 1);
        }
    }
    for(int i = 1; i < n; i++) sum[i] = add(sum[i - 1] + 1ll * i * i % p * phi[i] % p);
    inv6 = qpow(6, p - 2);
    return;
}

ll calc(ll n)
{
    if(n < N) return sum[n];
    if(mp[n]) return mp[n];
    ll res = s3(n);
    for(ll l = 2, r; l <= n; l = r + 1)
    {
        r = n / (n / l);
        res = sub(res - sub(s2(r) - s2(l - 1)) * calc(n / l) % p);
    }
    return mp[n] = res;
}

ll solve(ll n)
{
    ll res = 0;
    for(ll l = 1, r; l <= n; l = r + 1)
    {
        r = n / (n / l);
        res = add(res + sub(calc(r) - calc(l - 1)) * s3(n / l) % p);
    }
    return res;
}

int main()
{

    read(p), read(n);
    init(N);
    write(solve(n)), pc('\n');
    return 0;
}
// A.S.

Luogu P3704 数字表格

Description

TT 组数据，每组数据给定 n,mn,m，求

∏i=1n∏j=1mfgcd⁡(i,j)( mod  109+7)\prod_{i=1}^n\prod_{j=1}^mf_{\gcd(i,j)}\quad (\bmod\ 10^9+7)

其中 fif_i 表示斐波那契数列的第 ii 项。

1≤T≤103,1≤n,m≤1061\le T \le 10^3,1\le n,m\le 10^6

Solution

∏i=1n∏j=1mfgcd⁡(i,j)=∏i=1n∏j=1m∏d=1n[gcd⁡(i,j)=d]fd=∏d=1nfd∑i=1n∑j=1m[gcd⁡(i,j)=d]\begin{aligned} \prod_{i=1}^n\prod_{j=1}^mf_{\gcd(i,j)} &=\prod_{i=1}^n\prod_{j=1}^m\prod_{d=1}^n[\gcd(i,j)=d]f_{d} \\ &=\prod_{d=1}^nf_d^{\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=d]} \\ \end{aligned}

发现指数部分非常经典

∑i=1n∑j=1m[gcd⁡(i,j)=d]=∑k=1ndμ(k)⌊ndk⌋⌊mdk⌋\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)=d]=\sum_{k=1}^{\frac n d}\mu(k)\left\lfloor \frac n{dk} \right\rfloor\left\lfloor \frac m{dk} \right\rfloor

带回原式

∏d=1nfd∑k=1ndμ(k)⌊ndk⌋⌊mdk⌋=∏d=1n∏k=1ndfdμ(k)⌊ndk⌋⌊mdk⌋\prod_{d=1}^nf_d^{\sum_{k=1}^{\frac n d}\mu(k)\left\lfloor \frac n{dk} \right\rfloor\left\lfloor \frac m{dk} \right\rfloor} =\prod_{d=1}^n\prod_{k=1}^{\frac n d}f_d^{\mu(k)\left\lfloor \frac n{dk} \right\rfloor\left\lfloor \frac m{dk} \right\rfloor}

令 T=dkT=dk

∏T=1n∏d∣Tfdμ(Td)⌊ndk⌋⌊mdk⌋=∏T=1n(∏d∣Tfdμ(Td))⌊ndk⌋⌊mdk⌋\prod_{T=1}^n\prod_{d|T}f_d^{\mu(\frac T d)\left\lfloor \frac n{dk} \right\rfloor\left\lfloor \frac m{dk} \right\rfloor}=\prod_{T=1}^n\left(\prod_{d|T}f_d^{\mu(\frac T d)}\right)^{\left\lfloor \frac n{dk} \right\rfloor\left\lfloor \frac m{dk} \right\rfloor}

括号外面的可以整除分块，里面的可以直接预处理，复杂度是调和级数

Code 
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
105
#include <bits/stdc++.h>
#define ll long long
#define db double
#define gc getchar
#define pc putchar

using namespace std;

namespace IO
{
    template <typename T>
    void read(T &x)
    {
        x = 0; bool f = 0; char c = gc();
        while(!isdigit(c)) f |= c == '-', c = gc();
        while(isdigit(c)) x = x * 10 + c - '0', c = gc();
        if(f) x = -x;
    }

    template <typename T>
    void write(T x)
    {
        if(x < 0) pc('-'), x = -x;
        if(x > 9) write(x / 10);
        pc('0' + x % 10);
    }
}
using namespace IO;

const int N = 1e6 + 5;
const int mod = 1e9 + 7;

int p[N], tot, mu[N];
bool flag[N];
ll f[N], inv[N], pro[N];

ll add(ll x) {return x < mod ? x : x - mod;}

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

void init(int n)
{
    mu[1] = 1;
    f[1] = 1ll, inv[1] = 1ll, pro[0] = pro[1] = 1ll;

    for(int i = 2; i <= n; i++)
    {
        if(!flag[i]) p[++tot] = i, mu[i] = -1;
        for(int j = 1; j <= tot && i * p[j] <= n; j++)
        {
            flag[i * p[j]] = 1;
            if(i % p[j] == 0)
            {
                mu[i * p[j]] = 0;
                break;
            }
            mu[i * p[j]] = -mu[i];
        }
        f[i] = add(f[i - 1] + f[i - 2]);
        inv[i] = qpow(f[i], mod - 2);
        pro[i] = 1ll;
    }

    for(int i = 1; i <= n; i++)
        for(int j = i; j <= n; j += i)
           if(mu[j / i]) pro[j] = pro[j] * (mu[j / i] == 1 ? f[i] : inv[i]) % mod;
    for(int i = 2; i <= n; i++) pro[i] = pro[i] * pro[i - 1] % mod;

    return;
}

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

int main()
{
    init(1e6);
    int T; read(T);
    while(T--)
    {
        int n, m; read(n), read(m);
        if(n > m) n ^= m ^= n ^= m;
        write(solve(n, m)), pc('\n');
    }
    return 0;
}
// A.S.

Luogu P2257 YY的GCD

Description

TT 组数据，每组给定 n,mn,m，求

∑i=1n∑j=1m[gcd⁡(i,j)∈prime]\sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)\in prime]

T=104,n,m≤107T=10^4,n,m\le 10^7

Solution

默认 n≤mn\le m，P\mathbb{P} 表示质数

∑i=1n∑j=1m[gcd⁡(i,j)∈P]=∑i=1n∑j=1m∑d∈P[gcd⁡(i,j)=d]=∑d∈P∑i=1nd∑j=1md[gcd⁡(i,j)=1]=∑d∈P∑k=1ndμ(k)⌊ndk⌋⌊mdk⌋\begin{aligned} \sum_{i=1}^n\sum_{j=1}^m[\gcd(i,j)\in \mathbb{P}] &=\sum_{i=1}^n\sum_{j=1}^m\sum_{d\in \mathbb{P}}[\gcd(i,j)=d] \\ &=\sum_{d\in \mathbb{P}}\sum_{i=1}^{\frac n d}\sum_{j=1}^{\frac m d}[\gcd(i,j)=1] \\ &=\sum_{d\in \mathbb{P}}\sum_{k=1}^{\frac n d}\mu(k)\left\lfloor \frac n{dk} \right\rfloor\left\lfloor \frac m{dk} \right\rfloor \\ \end{aligned}

令 T=dkT=dk

∑T=1n∑d∣T,d∈Pμ(Td)⌊nT⌋⌊mT⌋=∑T=1n⌊nT⌋⌊mT⌋∑d∣T,d∈Pμ(Td)\sum_{T=1}^n\sum_{d|T,d\in \mathbb{P}}\mu(\frac T d)\left\lfloor \frac n{T} \right\rfloor\left\lfloor \frac m{T} \right\rfloor =\sum_{T=1}^n\left\lfloor \frac n{T} \right\rfloor\left\lfloor \frac m{T} \right\rfloor\sum_{d|T,d\in \mathbb{P}}\mu(\frac T d)

预处理一下后面的东西就可以整除分块了

Code 
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
#include <bits/stdc++.h>
#define ll long long
#define db double
#define gc getchar
#define pc putchar

using namespace std;

namespace IO
{
    template <typename T>
    void read(T &x)
    {
        x = 0; bool f = 0; char c = gc();
        while(!isdigit(c)) f |= c == '-', c = gc();
        while(isdigit(c)) x = x * 10 + c - '0', c = gc();
        if(f) x = -x;
    }

    template <typename T>
    void write(T x)
    {
        if(x < 0) pc('-'), x = -x;
        if(x > 9) write(x / 10);
        pc('0' + x % 10);
    }
}
using namespace IO;

const int N = 1e7 + 5;

int p[N], tot, mu[N], sum[N];
bool flag[N];

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

ll solve(int n, int m)
{
    if(n > m) n ^= m ^= n ^= m;
    ll res = 0;
    for(int l = 1, r; l <= n; l = r + 1)
    {
        r = min(n / (n / l), m / (m / l));
        res += 1ll * (sum[r] - sum[l - 1]) * (n / l) * (m / l);
    }
    return res;
}

int main()
{
    init(1e7);
    int T; read(T);
    while(T--)
    {
        int n, m;
        read(n), read(m);
        write(solve(n, m)), pc('\n');
    }
    return 0;
}
// A.S.