莫比乌斯反演例题集---(自用)

参考:大佬1 大佬2 大佬3 大佬4

P3455 [POI2007]ZAP-Queries

$$求解\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=k]$$ $反演过程:$ $$\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=k]$$

$$\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{m}{k} \right \rfloor}\varepsilon[gcd(i,j)=1]$$

$枚举d:$

$$\sum_{d=1}^{\left \lfloor min(\left \lfloor \frac{n}{k} \right \rfloor,\left \lfloor \frac{m}{k} \right \rfloor)\right \rfloor}\mu(d)\left \lfloor \frac{n}{kd} \right \rfloor \left \lfloor \frac{m}{kd} \right \rfloor$$

$\ {先预处理\mu,然后整数分块做,复杂度为O(n + T\sqrt n)}$

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
##include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

##define INF 0x7f7f7f
##define mem(a, b) memset(a , b , sizeof(a))
##define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 5e5 + 10;

int mu[N]; // 莫比乌斯函数
bool is\_prime[N];
int prime[N];
int cnt;

ll sum[N]; // 记录莫比乌斯函数前缀和

void Mobi() // 莫比乌斯函数初始化
{
mu[1] = 1;
is\_prime[0] = is\_prime[1] = true;
for(int i = 2;i < N; i++) {
if (!is\_prime[i]) {
mu[i] = -1;
prime[++cnt] = i;
}
for (int j = 1; j <= cnt && i \* prime[j] < N; j++) {
is\_prime[i \* prime[j]] = true;
if (i % prime[j] == 0) {
mu[i \* prime[j]] = 0;
break;
}
mu[i \* prime[j]] = -mu[i];
}
}

for(int i = 1;i < N; i++) {
sum[i] = sum[i - 1] + mu[i];
}
}

ll k;

ll Ans(ll n, ll m) {
ll ans = 0;
n /= k, m /= k;
int mx = min(n, m);
for(int l = 1, r;l <= mx; l = r + 1) {
r = min(n / (n / l), m / (m / l));
ans += (n / l) \* (m / l) \* (sum[r] - sum[l - 1]);
}
return ans;
}

void solve() {
Mobi();
int T;
cin >> T;
while(T--) {
ll n, m;
cin >> n >> m >> k;
cout << Ans(n, m) << endl;
}
}

signed main() {
ios\_base::sync\_with\_stdio(false);
//cin.tie(nullptr);
//cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

P2257 YY的GCD

$$求解\sum_{p\in prime}\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=p]$$ $反演过程和上面几乎一样,就多了一维p为质数的情况,直接化简为:$

$$\sum_{p\in prime}\sum_{d=1}^{\left \lfloor min(\left \lfloor \frac{n}{p} \right \rfloor,\left \lfloor \frac{m}{p} \right \rfloor)\right \rfloor}\mu(d)\left \lfloor \frac{n}{pd} \right \rfloor \left \lfloor \frac{m}{pd} \right \rfloor$$

$这里的复杂度为O(n+T\sqrt n \sqrt n),但明显不够,所以需要继续化简,这里的优化可以用T替换pd,则d=\frac{T}{p},替换后的式子为:$

$$\sum_{T=1}^{min(n,m)})\left \lfloor \frac{n}{T} \right \rfloor \left \lfloor \frac{m}{T} \right \rfloor\sum_{pT\;p\in prime}\mu(\frac{T}{p})$$

$$令\ {F(x)=\sum_{px\;p\in prime}\mu(\frac{x}{p})}$$

$得:$

$$\sum_{p\in prime}\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=p]=\sum_{T=1}^{min(n,m)}\left \lfloor \frac{n}{T} \right \rfloor \left \lfloor \frac{m}{T} \right \rfloor F(T)$$

$\ {O(n)处理\mu,O(nlogn)处理F,O(\sqrt n)询问分块,即复杂度为O(nlogn + T\sqrt n)}$

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
##include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

##define INF 0x7f7f7f
##define mem(a, b) memset(a , b , sizeof(a))
##define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 1e7 + 10;

int mu[N]; // 莫比乌斯函数
bool is\_prime[N];
int prime[N];
int cnt;

ll sum[N]; // 记录莫比乌斯函数前缀和
ll F[N];

void Mobi() // 莫比乌斯函数初始化
{
mu[1] = 1;
is\_prime[0] = is\_prime[1] = true;
for(int i = 2;i < N; i++) {
if (!is\_prime[i]) {
mu[i] = -1;
prime[++cnt] = i;
}
for (int j = 1; j <= cnt && i \* prime[j] < N; j++) {
is\_prime[i \* prime[j]] = true;
if (i % prime[j] == 0) {
mu[i \* prime[j]] = 0;
break;
}
mu[i \* prime[j]] = -mu[i];
}
}

for(int i = 1;i <= cnt; i++) {
for(int j = 1;j \* prime[i] < N; j++) {
F[j \* prime[i]] += mu[j];
}
}

for(int i = 1;i < N; i++) {
F[i] += F[i - 1];
}
}

ll k;

ll Ans(ll N, ll M) {
ll ans = 0;
ll n = N, m = M;
int mx = min(n, m);
for (int l = 1, r; l <= mx; l = r + 1) {
r = min(n / (n / l), m / (m / l));
ans += (n / l) \* (m / l) \* (F[r] - F[l - 1]);
}
return ans;
}

void solve() {
Mobi();
int T;
cin >> T;
while(T--) {
ll n, m;
cin >> n >> m;
cout << Ans(n, m) << endl;
}
}

signed main() {
ios\_base::sync\_with\_stdio(false);
//cin.tie(nullptr);
//cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

P1390 公约数的和

$$求解\sum_{i=1}^n\sum_{j=i+1}^ngcd(i,j)$$

$先来看看基本形式:$ $$\sum_{i=1}^n\sum_{j=1}^ngcd(i,j)$$ $则\ {ans=\frac{Calc(n,n)-\frac{n(n+1)}{2}}{2}},首先减去gcd(i,i)的个数\frac{n(n+1)}{2},然后减去gcd(i,j)=gcd(j,i)的个数,直接除2即可,则:$ $$\sum_{i=1}^n\sum_{j=1}^n\sum_{k=1}^nk[gcd(i,j)=k]$$

$枚举k:$

$$\sum_{k=1}^nk\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=\left \lfloor \frac{i+1}{k} \right \rfloor}^{\left \lfloor \frac{n}{k} \right \rfloor}[gcd(i,j)=1]$$

$$\sum_{k=1}^nk\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{di\;dj}\mu(d)$$

$$\sum_{k=1}^nk\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)\sum_{i=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}1$$

$$\sum_{k=1}^nk\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)(\left \lfloor \frac{n}{kd} \right \rfloor)^2$$

$令\ {T=kd并枚举T},得:$

$$\sum_{T=1}^n\sum_{dT}\mu(d)\frac{T}{d}(\left \lfloor \frac{n}{T} \right \rfloor)^2$$

$中间部分有\ {\sum_{dT}\mu(d)\frac{T}{d}=\varphi (T),因为\varphi =\mu * id},即得:$

$$\sum_{T=1}^n\varphi (T)(\left \lfloor \frac{n}{T} \right \rfloor)^2$$

$则最终答案为:$

$$\ {ans=\frac{\sum_{T=1}^n\varphi (T)(\left \lfloor \frac{n}{T} \right \rfloor)^2-\frac{n(n+1)}{2}}{2}}$$

$\ {O(n)预处理\varphi,\sqrt n询问分块,总时间复杂度为O(n+\sqrt n)}$

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
##include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

##define INF 0x7f7f7f
##define mem(a, b) memset(a , b , sizeof(a))
##define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const ll P = 1e9 + 7;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 2e6 + 10;

int phi[N];
bool is\_prime[N];
int prime[N];
int tot;

ll f[N];

void Euler()
{
phi[1] = 1;
is\_prime[1] = true;
for(int i = 2;i < N; i++){
if(!is\_prime[i]) {
phi[i] = i - 1;
prime[++tot] = i;
}
for(int j = 1;j <= tot && i \* prime[j] < N; j++){
is\_prime[i \* prime[j]] = true;
if(i % prime[j]) {
phi[i \* prime[j]] = phi[i] \* (prime[j] - 1);
}
else{
phi[i \* prime[j]] = phi[i] \* prime[j];
break;
}
}
}

for(int i = 1;i < N; i++) {
f[i] = f[i - 1] + phi[i];
}
}

ll Calc(ll n) {
ll ans = 0;

for(ll l = 1, r;l <= n;l = r + 1) {
r = min(n, n / (n / l));
ans += (f[r] - f[l - 1]) \* (n / l) \* (n / l);
}
return ans;
}

void solve() {
Euler();
ll n;
cin >> n;
cout << (Calc(n) - n \* (n + 1) / 2) / 2 << endl;
}

signed main() {
ios\_base::sync\_with\_stdio(false);
//cin.tie(nullptr);
//cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

P1829 [国家集训队]Crash的数字表格 / JZPTAB

$$求解\sum_{i=1}^n\sum_{j=1}^mlcm (i,j)$$

$由gcd和lcm之间的关系得:$

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

$$\sum_{k=1}^{min(n,m)}\sum_{i=1}^n\sum_{j=1}^m\frac{ij}{k}[gcd(i,j)=k]$$

$$\sum_{k=1}^{min(n,m)}\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{m}{k} \right \rfloor}\frac{ijk^2}{k}[gcd(i,j)=1]$$

$$\sum_{k=1}^{min(n,m)}k\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{m}{k} \right \rfloor}ij\sum_{di\;dj}\mu(d)$$

$枚举d得:$

$$\sum_{k=1}^{min(n,m)}k\sum_{d=1}^{min(\left \lfloor \frac{n}{k} \right \rfloor,\left \lfloor \frac{m}{k} \right \rfloor)}\mu(d)\sum_{i=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}id\sum_{j=1}^{\left \lfloor \frac{m}{kd} \right \rfloor}jd$$

$$\sum_{k=1}^{min(n,m)}k\sum_{d=1}^{min(\left \lfloor \frac{n}{k} \right \rfloor,\left \lfloor \frac{m}{k} \right \rfloor)}\mu(d)d^2\sum_{i=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}i\sum_{j=1}^{\left \lfloor \frac{m}{kd} \right \rfloor}j$$

$令\ {f(x)=\sum_{i=1}^{x}i,F(n,m)=\sum_{d=1}^{min(n,m)}\mu(d)d^2f(\frac{n}{d})f(\frac{m}{d})},得:$

$$\sum_{k=1}^{min(n,m)}kF(\left \lfloor\frac{n}{k}\right \rfloor,\left \lfloor\frac{m}{k}\right \rfloor)$$

$\ {O(n)处理\mu,O(1)处理f(等差数列,在计算F的过程中直接直接计算),O(\sqrt {\sqrt n})分块询问F,答案分块询问O(\sqrt n),总复杂度为O(n+n^{\frac{3}{4}})}$

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
##include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

##define INF 0x7f7f7f
##define mem(a, b) memset(a , b , sizeof(a))
##define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 1e7 + 10;

const ll mod = 20101009;

int mu[N]; // 莫比乌斯函数
bool is\_prime[N];
int prime[N];
int cnt;

ll sum[N];

void Mobi() // 莫比乌斯函数初始化
{
mu[1] = 1;
is\_prime[0] = is\_prime[1] = true;
for(int i = 2;i < N; i++) {
if (!is\_prime[i]) {
mu[i] = -1;
prime[++cnt] = i;
}
for (int j = 1; j <= cnt && i \* prime[j] < N; j++) {
is\_prime[i \* prime[j]] = true;
if (i % prime[j] == 0) {
mu[i \* prime[j]] = 0;
break;
}
mu[i \* prime[j]] = -mu[i];
}
}

for(int i = 1;i < N; i++) {
sum[i] = (sum[i - 1] + mu[i] \* i % mod \* i % mod) % mod;
}
}

ll Ans(ll n, ll m)
{
ll inv2 = (mod + 1) / 2;
ll ans = 0;
ll p = min(n, m);
for(int k = 1;k <= p; k++) {
int x = n / k, y = m / k, P = min(x, y);
ll Sum = 0;
for(int l = 1, r;l <= P; l = r + 1) {
r = min(x / (x / l), (y / (y / l)));
Sum = (Sum + (sum[r] - sum[l - 1] + mod) % mod \* (1 + x / l) % mod \* (x / l) % mod \* inv2 % mod \* (1 + y / l) % mod \* (y / l) % mod \* inv2 % mod) % mod;
}
ans = (ans + Sum \* k % mod) % mod;
}
return ans;
}

void solve() {
Mobi();
ll n, m;
cin >> n >> m;
cout << Ans(n, m) << endl;

}

signed main() {
ios\_base::sync\_with\_stdio(false);
//cin.tie(nullptr);
//cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

P3768 简单的数学题

$$求解\sum_{i=1}^n\sum_{j=1}^nij(gcd(i,j))\;mod \;p$$ $$\sum_{k=1}^n\sum_{i=1}^n\sum_{j=1}^nijk[gcd(i,j)=k]$$

$$\sum_{k=1}^nk\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{n}{k} \right \rfloor}ijk^2[gcd(i,j)=1]$$

$$\sum_{k=1}^nk^3\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{n}{k} \right \rfloor}ij\sum_{di\;dj}\mu(d)$$

$$\sum_{k=1}^nk^3\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)\sum_{i=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}id\sum_{j=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}jd$$

$$\sum_{k=1}^nk^3\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)d^2\sum_{i=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}i\sum_{j=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}j$$

$令\ {f(x)=\sum_{i=1}^xi},则式子变为:$

$$\sum_{k=1}^nk^3\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)d^2f^2(\frac{n}{kd})$$

$到这里,时间复杂度为O(n+n\sqrt n),我没尝试过,但是还可以继续化简,化简方法有替换,用T替换kd,则:$

$$\sum_{T=1}^nf^2(\frac{n}{T})T^2\sum_{kT}k*\mu(\frac{T}{k})$$

$$\sum_{T=1}^nf^2(\frac{n}{T})T^2\ {(id*\mu)T}$$

$$\sum_{T=1}^nf^2(\frac{n}{T})T^2\ {\varphi (T)}$$

$\ {可以杜教筛筛T^2\varphi (T)的前缀和,因为n最高有10^{10},明显数组存不下,所以对于5e6可以用数组直接过度求前缀和,对于后面就可以用杜教筛求前缀和,然后用map存下(STL的好处),过程如下:}$

$参考:$杜教筛

$先看杜教筛的基本形式(与上面无关):设S(n)=\sum_{i=1}^{n}f(i),定义一个新函数为g,则有$

$$\sum_{i=1}^n(f*g)i$$

$$\sum_{i=1}^n\sum_{xy=i}f(x)g(y)$$

$$\sum_{y=1}^ng(y)\sum_{xy\leq n}f(x)$$

$$\sum_{y=1}^ng(y)\sum_{x=1}^{\left \lfloor \frac{n}{y} \right \rfloor}f(x)$$

$$\sum_{y=1}^ng(y)S(\left \lfloor \frac{n}{y} \right \rfloor)$$

$则有:$

$$\sum_{i=1}^n(f*g)i=g(1)S(n)+\sum_{y=2}^ng(y)S(\left \lfloor \frac{n}{y} \right \rfloor)$$

$$\ {S(n)=\frac{\sum_{i=1}^n(f*g)i-\sum_{y=2}^ng(y)S(\left \lfloor \frac{n}{y} \right \rfloor)}{g(1)}}$$

$对于本题,要找到合适的g函数,因为我们求的前缀和为S(n)=\sum_{i=1}^ni^2\varphi(i),f(i)=i^2\varphi(i),所以设g(n)=n^2。$

$\sum_{i=1}^ng(i)=\frac{x(x+1)(2x+1)}{6}$

$(f*g)x=\sum_{dx}f(d)g(\frac{x}{d})=\sum_{dx}d^2\varphi(d)\frac{x^2}{d^2}=x^2\sum_{dx}\varphi(d)=x^2id(x)=x^3$

$\sum_{i=1}^n(f*g)i=\sum_{i=1}^ni^3=\frac{x^2(x+1)^2}{4}$

$$\ {S(n)=\frac{n^2(n+1)^2}{4}-\sum_{y=2}^ng(y)S(\left \lfloor \frac{n}{y} \right \rfloor)}$$

$\ {O(n^{\frac{2}{3}})处理T^2\varphi (T)的前缀和,f可以O(1)算出,O(\sqrt n) 询问分块,总时间复杂度为O(n^\frac{2}{3})}$

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
109
110
111
112
113
114
115
116
117
118
119
##include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

##define INF 0x7f7f7f
##define mem(a, b) memset(a , b , sizeof(a))
##define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 5e6 + 10;

ll n, mod, inv4, inv6;
ll phi[N];
bool is\_prime[N];
int prime[N];
int tot;

ll sum[N];

void Get\_phi() {
phi[1] = 1;
is\_prime[1] = true;
for(int i = 2;i <= N; i++){
if(!is\_prime[i]) {
phi[i] = i - 1;
prime[++tot] = i;
}
for(int j = 1;j <= tot && i \* prime[j] < N; j++){
is\_prime[i \* prime[j]] = true;
if(i % prime[j]) {
phi[i \* prime[j]] = phi[i] \* (prime[j] - 1);
}
else{
phi[i \* prime[j]] = phi[i] \* prime[j];
break;
}
}
}
for(int i = 1;i < N; i++) {
sum[i] = (sum[i - 1] + (ll)phi[i] \* i % mod \* i % mod) % mod;
}
}

ll quick\_pow(ll a, ll b) {
ll ans = 1;
while(b){
if(b & 1) ans = ans \* a % mod;
a = a \* a % mod;
b >>= 1;
}
return ans;
}

ll sum1(ll x) {x %= mod; return (1 + x) % mod \* x % mod \* (1 + x) % mod \* x % mod \* inv4 % mod;}
ll sum2(ll x) {x %= mod; return x % mod \* (1 + x) % mod \* (2 \* x + 1) % mod \* inv6 % mod;}

map<ll, ll> mp;

ll Calc(ll x) {
if(x < N) return sum[x];
if(mp[x]) return mp[x];
ll num = sum1(x);
for(ll l = 2, r;l <= x;l = r + 1) {
r = x / (x / l);
num = (num - Calc(x / l) % mod \* (sum2(r) - sum2(l - 1) + mod) % mod + mod) % mod;
}
return mp[x] = num;
}

ll Ans() {
ll ans = 0;
for(ll l = 1, r;l <= n; l = r + 1) {
r = n / (n / l);
ans = (ans + sum1(n / l) % mod \* (Calc(r) - Calc(l - 1) % mod + mod) % mod) % mod;
}
return ans;
}

void solve() {
cin >> mod >> n;
Get\_phi();
inv4 = quick\_pow(4, mod - 2);
inv6 = quick\_pow(6, mod - 2);
cout << Ans() << endl;
}

signed main() {
ios\_base::sync\_with\_stdio(false);
//cin.tie(nullptr);
//cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

P3327 [SDOI2015]约数个数和

$$求解\sum_{i=1}^n\sum_{j=1}^md (ij)$$

$$\ {约数的重要性质:(ij)=\sum_{xi}\sum_{yj}[gcd(x,y )=1]}$$

https://www.cnblogs.com/sun123zxy/p/12295533.html

$得:$ $$\sum_{i=1}^n\sum_{j=1}^m\sum_{xi}\sum_{yj}[gcd(x,y )=1]$$

$改变枚举顺序,枚举x和y:$

$$\sum_{x=1}^n\sum_{y=1}^m\left \lfloor \frac{n}{x} \right \rfloor\left \lfloor \frac{n}{y} \right \rfloor[gcd(x,y)=1]$$

$$\sum_{i=1}^n\sum_{j=1}^m\left \lfloor \frac{n}{i} \right \rfloor\left \lfloor \frac{n}{j} \right \rfloor[gcd(i,j)=1]$$

$后面套用前面得:$

$$\sum_{i=1}^n\sum_{j=1}^m\left \lfloor \frac{n}{i} \right \rfloor\left \lfloor \frac{n}{j} \right \rfloor\sum_{di\;dj}\mu(d)$$

$枚举d得:$

$$\sum_{d=1}^{min(n,m)}\mu(d)\sum_{i=1}^n\sum_{j=1}^m\left \lfloor \frac{n}{i} \right \rfloor\left \lfloor \frac{n}{j} \right \rfloor[dgcd(i,j)]$$

$由于[dgcd(i,j)]成立的条件为d是gcd(i,j)的约数,所以可以把枚举i,j变换为枚举di,dj,从而[1gcd(i,j)]可以省略。得:$

$$\sum_{d=1}^{min(n,m)}\mu(d)\sum_{i=1}^{\left \lfloor \frac{n}{d} \right \rfloor}\sum_{j=1}^{\left \lfloor\frac{m}{d}\right \rfloor}\left \lfloor \frac{n}{id} \right \rfloor\left \lfloor \frac{n}{jd} \right \rfloor$$

$$\sum_{d=1}^{min(n,m)}\mu(d)\sum_{i=1}^{\left \lfloor \frac{n}{d} \right \rfloor}\left \lfloor \frac{n}{id} \right \rfloor\sum_{j=1}^{\left \lfloor\frac{m}{d}\right \rfloor}\left \lfloor \frac{n}{jd} \right \rfloor$$

$令\ {f(x)=\sum_{i=1}^x\left \lfloor \frac{x}{i} \right \rfloor},则:$

$$\sum_{d=1}^{min(n,m)}\mu(d)f(\frac{n}{d})f(\frac{m}{d})$$

$\ {O(n)处理\mu和其前缀和,O(n\sqrt n)处理f,O(\sqrt n)询问分块,总复杂度为O(n\sqrt n+T\sqrt n)}$

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
##include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

##define INF 0x7f7f7f
##define mem(a, b) memset(a , b , sizeof(a))
##define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 5e4 + 10;

ll mu[N]; // 莫比乌斯函数
bool is\_prime[N];
int prime[N];
int cnt;

ll f[N];

void Mobi() // 莫比乌斯函数初始化
{
mu[1] = 1;
is\_prime[0] = is\_prime[1] = true;
for(int i = 2;i < N; i++) {
if (!is\_prime[i]) {
mu[i] = -1;
prime[++cnt] = i;
}
for (int j = 1; j <= cnt && i \* prime[j] < N; j++) {
is\_prime[i \* prime[j]] = true;
if (i % prime[j] == 0) {
mu[i \* prime[j]] = 0;
break;
}
mu[i \* prime[j]] = -mu[i];
}
}

for(int i = 1;i < N; i++) {
mu[i] += mu[i - 1];
}

for(int i = 1;i < N; i++) {
ll ans = 0;
for(int l = 1, r;l <= i; l = r + 1) {
r = i / (i / l);
ans += (r - l + 1) \* (i / l);
}
f[i] = ans;
}
}

ll Ans(ll n, ll m)
{
ll ans = 0;
ll k = min(n, m);

for(ll l = 1, r;l <= k; l = r + 1) {
r = min(n / (n / l), m / (m / l));
ans += (mu[r] - mu[l - 1]) \* f[m / l] \* f[n / l];
}
return ans;
}

void solve() {
int T;
Mobi();
cin >> T;
while(T--) {
ll n, m;
cin >> n >> m;
cout << Ans(n, m) << endl;
}
}

signed main() {
ios\_base::sync\_with\_stdio(false);
//cin.tie(nullptr);
//cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

P3312 [SDOI2014]数表

$$求解\sum_{i=1}^n\sum_{j=1}^m\sum_{di\;dj}d[\sum_{di\;dj}d\leq a]$$

$令\ {F(x)=\sum_{ix}i},则得:$

$$\sum_{i=1}^n\sum_{j=1}^m F(gcd(i,j))[F(gcd(i,j))\leq a]$$

$$\sum_{d=1}^{min(n,m)}\sum_{i=1}^n\sum_{j=1}^m F(d)[gcd(i,j)=d][F(d)\leq a]$$

$$\sum_{d=1}^{min(n,m)}F(d)[F(d)\leq a]\sum_{i=1}^n\sum_{j=1}^m[gcd(i,j)=d]$$

$后半部分就套用前面得:$

$$\sum_{d=1}^{min(n,m)}F(d)[F(d)\leq a]\sum_{d^{‘}=1}^{min(\left \lfloor \frac{n}{d} \right \rfloor, \left \lfloor \frac{m}{d} \right \rfloor)}\mu(d^{‘})\left \lfloor \frac{n}{dd^{‘}} \right \rfloor \left \lfloor \frac{m}{dd^{‘}} \right \rfloor$$

$令T替换dd^{‘}得:$

$$\sum_{T=1}^{min(n,m)}\sum_{dT}F(d)\mu(\frac{T}{d})\left \lfloor \frac{n}{T} \right \rfloor \left \lfloor \frac{m}{T} \right \rfloor[F(d)\leq a]$$

$$\sum_{T=1}^{min(n,m)}\left \lfloor \frac{n}{T} \right \rfloor \left \lfloor \frac{m}{T} \right \rfloor\sum_{dT}F(d)\mu(\frac{T}{d})[F(d)\leq a]$$

$令\ {f(x)=\sum_{dx}F(d)\mu(\frac{x}{d})[F(d)\leq a]}得:$

$$\sum_{T=1}^{min(n,m)}\left \lfloor \frac{n}{T} \right \rfloor \left \lfloor \frac{m}{T} \right \rfloor f(T)$$

$\ {O(n)处理\mu,O(nlogn)处理F,对于每次询问的a值可以离线排序,树状数组维护f,处理时间为O(nlog^{2}n),O(\sqrt n\; logn)询问分块,总时间复杂度为O(nlog^2n)+T\sqrt n\; logn}$

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
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
##include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

##define INF 0x3f3f3f3f
##define lc u << 1
##define rc u << 1 1
##define m (l + r) / 2
##define mid (t[u].l + t[u].r) / 2
##define lowbit(x) x & (-x)
##define mem(a, b) memset(a , b , sizeof(a))
##define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const ll mod = 1e9 + 7;
// const double eps = 1e-6;
// const double PI = acos(-1);
// const double R = 0.57721566490153286060651209;

const int N = 1e5 + 1;
ll mod = (1ll \* 1) << 31;

struct node {
int id;
ll data;
}F[N];
struct Node {
int nn, mm, id;
ll a;
}q[N];

int cnt, prime[N];
bool is\_prime[N];
int mu[N];

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

bool cmp\_Q(Node a, Node b) {
return a.a < b.a;
}

inline void init()
{
mu[1] = 1;
is\_prime[0] = is\_prime[1] = true;
for(int i = 2;i < N; i++) {
if (!is\_prime[i]) {
mu[i] = -1;
prime[++cnt] = i;
}
for (int j = 1; j <= cnt && i \* prime[j] < N; j++) {
is\_prime[i \* prime[j]] = true;
if (i % prime[j] == 0) {
mu[i \* prime[j]] = 0;
break;
}
mu[i \* prime[j]] = -mu[i];
}
}

for(int i = 1;i < N; i++) {
F[i].id = i;
for(int j = 1;j \* i < N; j++) {
F[i \* j].data += i;
}
}
sort(F + 1, F + N, cmp\_F);

}

ll t[N];
ll ans[N];

inline void add(int x, int v) {
while(x < N) {
t[x] = (t[x] + v) % mod;
if(t[x] > mod)
t[x] %= mod;
x += lowbit(x);
}
}

inline ll query(int x) {
ll ans = 0;
while(x) {
ans = (ans + t[x]) % mod;
if(ans > mod)
ans %= mod;
x -= lowbit(x);
}
return ans;
}

inline ll Calc(int nn, int mm) {
if(nn > mm)
swap(nn, mm);

ll ans = 0;
for(int l = 1, r;l <= nn; l = r + 1) {
r = min(nn / (nn / l), mm / (mm / l));
ans = (ans + 1ll \* (nn / l) \* (mm / l) % mod \* (query(r) - query(l - 1)) % mod) % mod;
}
return ans;
}

//struct SegMent {
// int l, r;
// int sum;
//}t[N << 2];
//
//void push\_up(int u) {
// t[u].sum = t[lc].sum + t[rc].sum;
//}
//
//void build(int u, int l, int r) {
// t[u].l = l;
// t[u].r = r;
// if(l == r) {
// t[u].sum = a[l];
// return ;
// }
// build(lc, l, m);
// build(rc, m + 1, r);
// push\_up(u);
//}
//
//void update(int u, int p, int v) {
// if(t[u].l == t[u].r) {
// t[u].sum += v;
// return ;
// }
// if(p <= mid) update(lc, p, v);
// else update(rc, p, v);
// push\_up(u);
//}
//
//int query(int u, int ql, int qr) {
// if(ql <= t[u].l && t[u].r <= qr) {
// return t[u].sum;
// }
// int ans = 0;
// if(ql <= mid)
// ans += query(lc, ql, qr);
// if(qr > mid)
// ans += query(rc, ql, qr);
// return ans;
//}

void solve() {
init();
int T;
scanf("%d",&T);
for(int i = 1;i <= T; i++) {
scanf("%d%d%lld",&q[i].nn,&q[i].mm, &q[i].a);
q[i].id = i;
}
sort(q + 1, q + T + 1, cmp\_Q);

int now = 1;
for(int i = 1;i <= T; i++) {
while(now < N && F[now].data <= q[i].a) {
for(int j = 1;j \* F[now].id < N; j++) {
add(j \* F[now].id, mu[j] \* F[now].data);
}
now++;
}
ans[q[i].id] = (Calc(q[i].nn, q[i].mm) % mod + mod) % mod;
}
for(int i = 1;i <= T; i++) {
printf("%lld\n",ans[i]);
}
}

signed main() {
ios\_base::sync\_with\_stdio(false);
//cin.tie(nullptr);
//cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

P3172 [CQOI2015]选数

$$求解\sum_{i=L}^H\sum_{j=L}^H…\sum_{x=L}^H[gcd(i,j…x)=k]$$

$本题求n个数的gcd为k的个数,根据上面套路,将H/k,L/k,这样就只用求n个数的gcd为1的个数了,注意:如果L整除k,则L=L/k,否则L=L/k+1,就是向上取整,那么有一种办法非常好,就让\ {L/k+(\frac{k-1}{k})},即:$

$$\sum_{i={\left \lfloor \frac{L-1}{k} \right \rfloor+1}}^{\left \lfloor \frac{H}{k}\right \rfloor}\sum_{j={\left \lfloor \frac{L-1}{k} \right \rfloor+1}}^{\left \lfloor \frac{H}{k}\right \rfloor}…\sum_{x={\left \lfloor \frac{L-1}{k} \right \rfloor+1}}^{\left \lfloor \frac{H}{k}\right \rfloor}[gcd(i,j…x)=1]$$

$$\sum_{i={\left \lfloor \frac{L-1}{k} \right \rfloor+1}}^{\left \lfloor \frac{H}{k}\right \rfloor}\sum_{j={\left \lfloor \frac{L-1}{k} \right \rfloor+1}}^{\left \lfloor \frac{H}{k}\right \rfloor}…\sum_{x={\left \lfloor \frac{L-1}{k} \right \rfloor+1}}^{\left \lfloor \frac{H}{k}\right \rfloor}\sum_{di\;dj\;..\;dx}\mu(d)$$

$枚举d,令\ {l=\left \lfloor \frac{L-1}{k} \right \rfloor+1,r=\left \lfloor \frac{H}{k}\right \rfloor}:$

$$\sum_{d=1}^r\mu(d)\sum_{i=l}^r[di]\sum_{j=l}^r[dj]…\sum_{x=l}^r[dx]$$

$只有当i,j…x都是d的倍数的时候,最终答案+1,根据容斥原理,那么在[l,r]中d倍数有\ {\frac{r}{d}-\frac{l-1}{d}},对于所有的i,j…x,他们当中所有[l,r]的d的倍数的个数有\ {(\frac{r}{d}-\frac{l-1}{d})^n},即为:$

$$\sum_{d=1}^{r}\mu(d)*(\frac{r}{d}-\frac{l-1}{d})^n$$

$这里的H非常大,如果用线性筛去处理\mu,会直接Wa掉,因为没法处理到10^9次方那么大的前缀和,数组也开不了那么大。后果:$

在这里插入图片描述 $这里的AC可能就是一小部分的运气罢了。如果用线性筛,RE的RE,WA的WA,所以推荐用杜教筛取筛。$

$\ {设f(d)=\mu(d),S(r)=\sum_{d=1}^nf(d),这里的S就是要求的前缀和。}$ $我们选择\ {g=I},则\sum_{i=1}^r(f*g)i=\sum_{i=1}^r(\mu*I)i=\sum_{i=1}^r\varepsilon (i)=1。$

$枚举y:$

$$\sum_{y=1}^rg(y)\sum_{x=1}^{\left \lfloor \frac{r}{y} \right \rfloor)}f(x)$$

$$\sum_{y=1}^rg(y)S(\left \lfloor \frac{r}{y} \right \rfloor)$$

$$\sum_{i=1}^r(f*g)i =g(1)S(r)+\sum_{y=2}^rg(y)S(\left \lfloor \frac{r}{y} \right \rfloor)$$

$$S(r)=\frac{\sum_{i=1}^r(f*g)i-\sum_{y=2}^rg(y)S(\left \lfloor \frac{r}{y} \right \rfloor)}{g(1)}$$

$当g=I时,即:$

$$\ {S(r)=1-\sum_{y=2}^rS(\left \lfloor \frac{r}{y} \right \rfloor)}$$

$\ {杜教筛的复杂度为O(n\frac{2}{3}),比线性筛还要快一点,总时间复杂度为O(n^{\frac{2}{3}})。}$

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
109
110
111
112
113
114
115
116
117
118
119
##include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

##define INF 0x7f7f7f
##define mem(a, b) memset(a , b , sizeof(a))
##define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const ll mod = 1e9 + 7;

const int N = 1e6 + 10;

int mu[N]; // 莫比乌斯函数
bool is\_prime[N];
int prime[N];
int cnt;
ll sum[N];

ll n, L, H, k;

ll quick\_pow(ll a, ll b) {
ll ans = 1;
while(b) {
if(b & 1) ans = ans \* a % mod;
a = a \* a % mod;
b >>= 1;
}
return ans;
}

void Mobi() // 莫比乌斯函数初始化
{
mu[1] = 1;
is\_prime[0] = is\_prime[1] = true;
for(int i = 2;i < N; i++) {
if (!is\_prime[i]) {
mu[i] = -1;
prime[++cnt] = i;
}
for (int j = 1; j <= cnt && i \* prime[j] < N; j++) {
is\_prime[i \* prime[j]] = true;
if (i % prime[j] == 0) {
mu[i \* prime[j]] = 0;
break;
}
mu[i \* prime[j]] = -mu[i];
}
}
for(int i = 1;i < N; i++) {
sum[i] = (sum[i - 1] + mu[i]) % mod;
}
}

map<ll , ll> mp;

ll Calc(ll x) {
if(x < N)
return sum[x];
if(mp[x])
return mp[x];
ll num = 1;
for(int l = 2, r;l <= x; l = r + 1) {
r = x / (x / l);
num = (num - Calc(x / l) \* (r - l + 1) % mod + mod) % mod;
}
return mp[x] = num;
}

ll Ans() {
ll ans = 0;
for(ll l = 1, r;l <= H;l = r + 1) {
r = min(H / (H / l), (L / l) ? L / (L / l) : H + 2);
ans = (ans + (Calc(r) - Calc(l - 1) + mod) % mod \* quick\_pow((H / l) - (L / l), n) % mod + mod) % mod;
}
return (ans + mod) % mod;
}

void solve() {
cin >> n >> k >> L >> H;
L = (L - 1) / k; // 因为后面需要L-1,所以这里干脆就不+1了,但是含义不同
H = H / k;
Mobi();
cout << Ans() << endl;
}

signed main() {
ios\_base::sync\_with\_stdio(false);
//cin.tie(nullptr);
//cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

AT5200 [AGC038C] LCMs

$$求解\sum_{i=1}^{N}\sum_{j=1}^Nlcm(A_i,A_j)$$

$莫比乌斯反演处理的一般都是枚举变量之间的关系,所以我们要让A_i和A_j变为i和j的形式。$

$想要用i和j表示A_i和A_j并且不遗漏任何的A,所以让A等于某个i,就是让cnt_i=\ {\sum_{d=1}^N[A_d=i]},统计所有A等于某个i的个数,即cnt_i,j也同样如此,如下所示:$

$$\sum_{i=1}^{N}\sum_{j=1}^Ncnt_icnt_jlcm(i,j)$$

$这样还不够,因为i只能是[1,N]的,当Ai中有些值大于N的话,再枚举i从[1,N]就会漏掉大于N的Ai的情况,所以想要统计所有情况,即让N=Ai中最大的数即可,即\ {M=\max_{i=1}^N }。$ $所以我们需要反演的式子为:$

$$\sum_{i=1}^{M}\sum_{j=1}^Mcnt_i ⋅cnt_j ⋅lcm(i,j)$$

$$\sum_{i=1}^{M}\sum_{j=1}^Mcnt_i ⋅cnt_j ⋅\frac{ij}{gcd(i,j)}$$

$$\sum_{k=1}^M\sum_{i=1}^{M}\sum_{j=1}^Mcnt_i ⋅cnt_j ⋅\frac{ij}{k[gcd(i,j)=k]}$$

$$\sum_{k=1}^M\sum_{i=1}^{\left \lfloor \frac{M}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{M}{k} \right \rfloor}cnt_{ik} ⋅cnt_{jk} ⋅\frac{ijk^2}{k}[gcd(i,j)=1]$$

$$\sum_{k=1}^Mk\sum_{i=1}^{\left \lfloor \frac{M}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{M}{k} \right \rfloor}cnt_{ik} ⋅cnt_{jk} ⋅i ⋅j\sum_{di\;dj}\mu(d)$$

$$\sum_{k=1}^Mk\sum_{d=1}^{\left \lfloor \frac{M}{k} \right \rfloor}\mu(d)\sum_{i=1}^{\left \lfloor \frac{M}{kd} \right \rfloor}cnt_{ikd} ⋅id\sum_{j=1}^{\left \lfloor \frac{M}{kd} \right \rfloor}cnt_{jkd} ⋅jd$$

$$\sum_{k=1}^Mk\sum_{d=1}^{\left \lfloor \frac{M}{k} \right \rfloor}\mu(d)d^2(\sum_{i=1}^{\left \lfloor \frac{M}{kd} \right \rfloor}cnt_{ikd} ⋅i)^2$$

$令T=kd则:$

$$\sum_{T=1}^M(\sum_{i=1}^{\left \lfloor \frac{M}{T} \right \rfloor}cnt_{iT} ⋅i)^2\sum_{kT}\mu(k)k^2\frac{T}{k}$$

$$\sum_{T=1}^MT(\sum_{i=1}^{\left \lfloor \frac{M}{T} \right \rfloor}cnt_{iT} ⋅i)^2\sum_{kT}\mu(k)k$$

$令\ {f(x)=\sum_{ix}\mu(i)i\;,\;g(x)=\sum_{i=1}^{\left \lfloor \frac{M}{x} \right \rfloor}cnt_{ix} ⋅i=\frac{1}{x}\sum_{xt}tcnt_t},则最后的式子为:$

$$\sum_{T=1}^MTg^2(T)f(T)$$

$\ {O(n)处理A的次数,O(nlogn)处理f和g,总时间复杂度为O(nlogn)}$

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>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

##define INF 0x7f7f7f
##define mem(a, b) memset(a , b , sizeof(a))
##define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 5e5 + 10;

int mu[N]; // 莫比乌斯函数
bool is\_prime[N];
int prime[N];
int tot;
ll f[N], g[N];

int M, n;
int cnt[N];

void Mobi() // 莫比乌斯函数初始化
{
mu[1] = 1;
is\_prime[0] = is\_prime[1] = true;
for(int i = 2;i < N; i++) {
if (!is\_prime[i]) {
mu[i] = -1;
prime[++tot] = i;
}
for (int j = 1; j <= tot && i \* prime[j] < N; j++) {
is\_prime[i \* prime[j]] = true;
if (i % prime[j] == 0) {
mu[i \* prime[j]] = 0;
break;
}
mu[i \* prime[j]] = -mu[i];
}
}
for(int i = 1;i <= M; i++) {
for(int x = i;x <= M; x += i) {
f[x] += mu[i] \* i;
}
}

for(int x = 1;x <= M; x++) {
for(int t = x;t <= M; t += x) {
g[x] += t \* cnt[t];
}
g[x] /= x;
}
}

void solve() {
cin >> n;
M = n;
for(int i = 1;i <= n; i++) {
int x;
cin >> x;
cnt[x]++;
M = max(M, x);
}
Mobi();

ll ans = 0;

for(int i = 1;i <= M; i++) {
ans += i \* g[i] \* g[i] \* f[i];
}

cout << ans << endl;
}

signed main() {
ios\_base::sync\_with\_stdio(false);
//cin.tie(nullptr);
//cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

SP5971 LCMSUM - LCM Sum

$$求解\sum_{i=1}^nlcm(i,n)$$

$\left \lfloor \frac{n}{k} \right \rfloor$

$$\sum_{i=1}^n\frac{in}{gcd(i,n)}$$

$$\sum_{k=1}^n\sum_{i=1}^n\frac{in}{k[gcd(i,n)=k]}$$

$$\sum_{k=1}^n\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}in[gcd(i,\frac{n}{k})=1]$$

$$n\sum_{k=1}^n\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}i[gcd(i,\frac{n}{k})=1]$$

$$n\sum_{kn}\sum_{i=1}^ki[gcd(i,k)=1]$$

$$n\sum_{kn}\sum_{i=1}^ki\sum_{di\;dk}\mu(d)$$

$枚举d:$

$$n\sum_{kn}\sum_{dk}\mu(d)\sum_{i=1}^{\left \lfloor \frac{k}{d} \right \rfloor}id$$

$$n\sum_{kn}\sum_{dk}\mu(d)d\sum_{i=1}^{\left \lfloor \frac{k}{d} \right \rfloor}i$$

$$n\sum_{kn}\sum_{dk}\mu(d)d\frac{\left \lfloor \frac{k}{d} \right \rfloor*(\left \lfloor \frac{k}{d} \right \rfloor+1)}{2}$$

$$\frac{n}{2}\sum_{kn}{\sum_{dk}\mu(d)*\frac{k^2}{d}+\sum_{dk}\mu(d)*k}$$

$$\frac{n}{2}\sum_{kn}k{\sum_{dk}\mu(d)*\frac{k}{d}+\sum_{dk}\mu(d)}$$

$\ {这时候发现两个狄利克雷卷积,\varphi (n)=\mu*id=\sum_{dn}\mu(d)*\frac{n}{d}\;\;,\;\;\varepsilon(n) =\mu*I=\sum_{dn}\mu(d)。所以我们将其替换得:}$

$$\frac{n}{2}\sum_{kn}k[\varphi(k)+\varepsilon(k)]$$

$\ {O(n)处理\varphi,\sqrt n询问答案,总复杂度为O(n+T\sqrt n)。}$

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
##include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

##define INF 0x7f7f7f
##define mem(a, b) memset(a , b , sizeof(a))
##define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

const int N = 1e6 + 10;

int phi[N];
bool is\_prime[N];
int prime[N];
int tot;

ll f[N];

void Euler()
{
phi[1] = 1;
is\_prime[1] = true;
for(int i = 2;i < N; i++){
if(!is\_prime[i]) {
phi[i] = i - 1;
prime[++tot] = i;
}
for(int j = 1;j <= tot && i \* prime[j] < N; j++){
is\_prime[i \* prime[j]] = true;
if(i % prime[j]) {
phi[i \* prime[j]] = phi[i] \* (prime[j] - 1);
}
else{
phi[i \* prime[j]] = phi[i] \* prime[j];
break;
}
}
}
for(int i = 1;i < N; i++) {
for(int j = i;j < N; j += i) {
f[j] += 1ll \* phi[i] \* i + (i == 1);
}
}
}

void solve() {
Euler();
int T;
cin >> T;
while(T--) {
int n;
cin >> n;
cout << f[n] \* n / 2 << endl;
}
}

signed main() {
ios\_base::sync\_with\_stdio(false);
//cin.tie(nullptr);
//cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

P5221 Product

$$求解\prod_{i=1}^n\prod_{j=1}^n\frac{lcm(i,j)}{gcd(i,j)}$$

$首先还是把lcm(i,j)化简得:$

$$\prod_{i=1}^n\prod_{j=1}^n\frac{ij}{[gcd(i,j)]^2}$$

$$\prod_{i=1}^n\prod_{j=1}^nij*\prod_{i=1}^n\prod_{j=1}^n\frac{1}{[gcd(i,j)]^2}$$

$\ {我们很容易推出前半部分的数值为(n!)^{2n },后半部分的值就是\prod_{i=1}^n\prod_{j=1}^n[gcd(i,j)]^2的逆元的平方。}$

$所以现在求解:$

$$\prod_{i=1}^n\prod_{j=1}^ngcd(i,j)$$

$$\prod_{k=1}^n\prod_{i=1}^n\prod_{j=1}^nk[gcd(i,j)=k]$$

$$\prod_{k=1}^nk^{\sum_{i=1}^n\sum_{j=1}^n[gcd(i,j)=k]}$$

$这个式子的变换,本来是[1,n]中所有的gcd累乘,然后我们单独把gcd(i,j)=k的拿出来,k会有[1,n]的范围,然后我们再将这些k累乘就会达到之前的效果。$ $不过,[1,n]中的gcd不单单只有一个,可能会有多个相同的gcd,所以变换之后的式子里还需要统计gcd=k出现的次数,那只需要将他们累加即可。$

$所以没毛病。继续:$

$$\prod_{k=1}^nk^{\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{n}{k} \right \rfloor}[gcd(i,j)=1]}$$

$$\prod_{k=1}^nk^{\sum_{i=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\sum_{di\;dj}\mu(d)}$$

$枚举d:$

$$\prod_{k=1}^nk^{\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)\sum_{i=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}\sum_{j=1}^{\left \lfloor \frac{n}{kd} \right \rfloor}1}$$

$$\ {\prod_{k=1}^nk^{\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)(\left \lfloor \frac{n}{kd} \right \rfloor)^2}}$$

$最终答案为:$

$$\ {(n!)^{2n}\prod_{k=1}^nk^{\sum_{d=1}^{\left \lfloor \frac{n}{k} \right \rfloor}\mu(d)(\left \lfloor \frac{n}{kd} \right \rfloor)^2}}$$

$\ {O(n)处理\mu的前缀和,\sqrt n询问指数的分块,可以用欧拉降幂优化,总时间复杂度为O(n+nlog_2n}$

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>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

##define INF 0x7f7f7f
##define mem(a, b) memset(a , b , sizeof(a))
##define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
// const ll P = 1e9 + 7;
// const int maxn = 1e5 + 10;
// const double eps = 1e-6;

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

bool is\_prime[N];
int prime[N];
int mu[N]; // 莫比乌斯函数
int cnt;

int n;
ll fac = 1;

ll quick\_pow(ll a, ll b)
{
ll ans = 1;
while(b) {
if(b & 1) ans = ans \* a % mod;
a = a \* a % mod;
b >>= 1;
}
return ans;
}
void Mobi() // 莫比乌斯函数初始化
{
mu[1] = 1;
is\_prime[0] = is\_prime[1] = true;
for(int i = 2;i <= n; i++) {
if (!is\_prime[i]) {
mu[i] = -1;
prime[++cnt] = i;
}
for (int j = 1; j <= cnt && i \* prime[j] <= n; j++) {
is\_prime[i \* prime[j]] = true;
if (i % prime[j] == 0) {
mu[i \* prime[j]] = 0;
break;
}
mu[i \* prime[j]] = -mu[i];
}
mu[i] = mu[i - 1] + mu[i];
fac = fac \* i % mod;
}
}

ll Ans() {
ll ans = 1;
for(int k = 1;k <= n; k++) {
ll p = 0;
int t = n / k;
for(int l = 1, r;l <= t; l = r + 1) {
r = min(t, t / (t / l));
p = (p + (mu[r] - mu[l - 1] + mod - 1) \* (t / l) % (mod - 1) \* (t / l) % (mod - 1)) % (mod - 1);
}
ans = ans \* quick\_pow(k, p) % mod;
}
ll k = quick\_pow(ans, mod - 2);
fac = quick\_pow(fac, n);
return quick\_pow(k, 2) \* fac % mod \* fac % mod;
}

void solve() {
cin >> n;
Mobi();
cout << Ans() << endl;
}

signed main() {
ios\_base::sync\_with\_stdio(false);
//cin.tie(nullptr);
//cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

HDU 6588 Function

$$求解\sum_{i=1}^ngcd(\left \lfloor \sqrt[3]i \right \rfloor,i)$$

$我们先将这个一重和式转化为二重和式,得:$

$$\sum_{a=1}^{\left \lfloor \sqrt[3]n \right \rfloor}\sum_{i=1}^ngcd(a,i)[\left \lfloor \sqrt[3]i \right \rfloor=a]$$

$我们来看这个式子,只有当\left \lfloor \sqrt[3]i \right \rfloor=a时,才会有对答案的贡献为gcd(a,i),所以\ {a \leq \left \lfloor \sqrt[3]i \right \rfloor < a+1,a^3 \leq i \leq (a+1)^3-1},所以上式转化为:$

$$\ {\sum_{a=1}^{\left \lfloor \sqrt[3]n \right \rfloor}\sum_{i=a^3}^{min(n,(a+1)^3-1)}gcd(a,i)}$$

$然后我们将这个二重和式拆开得:$

$${\color{Brown}\sum_{a=1}^{\left \lfloor \sqrt[3]n \right \rfloor-1}\sum_{i=a^3}^{min(n,(a+1)^3-1)}gcd(a,i)+\sum_{i=(\left \lfloor \sqrt[3]n \right \rfloor)^3}^ngcd(\left \lfloor \sqrt[3]n \right \rfloor, i)}$$

$令\sqrt[3]n=N,并使a=i,i=j(好看一点)得:$

$${\color{Brown}\sum_{i=1}^{N-1}\sum_{j=i^3}^{(i+1)^3-1}gcd(i,j)+\sum_{j=N^3}^ngcd(N, j)}$$

$先计算前半部分:$

$——————————–$

HDU 6833 A Very Easy Math Problem

$$求解\sum_{a_1=1}^n\sum_{a_2=1}^n…\sum_{a_x=1}^n\left(\prod_{j=1}^xa_j^k\right)f(gcd(a_1,a_2…a_x))gcd(a_1,a_2…a_x)$$ $令d=gcd(a_1,a_2…a_x)得:$

$$\sum_{a_1=1}^n\sum_{a_2=1}^n…\sum_{a_x=1}^n\left(\prod_{j=1}^xa_j^k\right)f(d)d$$

$枚举d得:$

$$\sum_{d=1}^ndf(d)\sum_{a_1=1}^n\sum_{a_2=1}^n…\sum_{a_x=1}^n\left(\prod_{j=1}^xa_j^k\right)[gcd(a_1,a_2…a_x)=d]$$

$$\sum_{d=1}^nd^{kx+1}f(d)\sum_{a_1=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\sum_{a_2=1}^{\left \lfloor \frac{n}{d} \right \rfloor }…\sum_{a_x=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\left(\prod_{j=1}^xa_j^k\right)[gcd(a_1,a_2…a_x)=1]$$

$$\sum_{d=1}^nd^{kx+1}f(d)\sum_{a_1=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\sum_{a_2=1}^{\left \lfloor \frac{n}{d} \right \rfloor }…\sum_{a_x=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\left(\prod_{j=1}^xa_j^k\right)\sum_{ta_1\;ta2\;…\;ta_x}\mu(t)$$

$$\sum_{d=1}^nd^{kx+1}f(d)\left (\sum_{i=1}^{\left \lfloor \frac{n}{d} \right \rfloor }i^k\right)^x\sum_{ta_1\;ta2\;…\;ta_x}\mu(t)$$

$枚举t得:$

$$\sum_{d=1}^nd^{kx+1}f(d)\sum_{t=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\mu(t)\left (\sum_{i=1}^{\left \lfloor \frac{n}{dt} \right \rfloor }(it)^k\right)^x$$

$$\sum_{d=1}^nd^{kx+1}f(d)\sum_{t=1}^{\left \lfloor \frac{n}{d} \right \rfloor }\mu(t)t^{kx}\left (\sum_{i=1}^{\left \lfloor \frac{n}{dt} \right \rfloor }i^k\right)^x$$

$令T=dt并枚举T得:$

$$\sum_{T=1}^n\left (\sum_{i=1}^{\left \lfloor \frac{n}{dt} \right \rfloor }i^k\right)^xT^{kx}\sum_{dT}df(d)\mu(\frac{T}{d})$$

$令g(T)=\sum_{dT}df(d)\mu(\frac{T}{d}),则$

$$\sum_{T=1}^n\left (\sum_{i=1}^{\left \lfloor \frac{n}{T} \right \rfloor }i^k\right)^xT^{kx}g(T)$$

$\ {O(nlogn)筛f和g,并对T^{kx}g(T)和\sum_{i=1}^{\left \lfloor \frac{n}{T} \right \rfloor }i^k做前缀和,最后\sqrt n分块计算,复杂度为O(nlogn+T\sqrt n)}$

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
109
110
111
112
113
114
115
116
117
118
119
120
121

##include <bits/stdc++.h>

using namespace std;

typedef long long ll;
typedef long double ld;
typedef pair<int, int> pdd;

##define INF 0x3f3f3f3f
##define lowbit(x) x & (-x)
##define mem(a, b) memset(a , b , sizeof(a))
##define FOR(i, x, n) for(int i = x;i <= n; i++)

// const ll mod = 998244353;
const ll mod = 1e9 + 7;
// const double eps = 1e-6;
// const double pi = acos(-1);

const int N = 2e5 + 10;

ll T, k, x;

int mu[N]; // 莫比乌斯函数
bool is\_prime[N];
int prime[N];
int cnt;
ll g[N], f[N];
ll sumG[N], sumi[N];

ll quick\_pow(ll a, ll b) {
ll ans = 1;
while(b) {
if(b & 1) ans = ans \* a % mod;
a = a \* a % mod;
b >>= 1;
}
return ans % mod;
}

void Mobi() // 莫比乌斯函数初始化
{
f[1] = mu[1] = 1;
is\_prime[0] = is\_prime[1] = true;
for(int i = 2;i < N; i++) {
f[i] = 1;
if (!is\_prime[i]) {
mu[i] = -1;
prime[++cnt] = i;
}
for (int j = 1; j <= cnt && i \* prime[j] < N; j++) {
is\_prime[i \* prime[j]] = true;
if (i % prime[j] == 0) {
mu[i \* prime[j]] = 0;
break;
}
mu[i \* prime[j]] = -mu[i];
}
}
// 筛f函数
for(int d = 2;d \* d < N; d++) {
for(int i = d \* d;i < N; i += d \* d) {
f[i] = 0;
}
}
// 筛g函数
for(int d = 1;d < N; d++) {
for(int i = d;i < N; i += d) {
g[i] = (g[i] + 1ll \* d \* f[d] % mod \* mu[i / d] % mod + mod) % mod;
}
}
// 处理sumG函数和sumi函数前缀和
for(int i = 1;i < N; i++) {
ll t = quick\_pow(i, k);
sumi[i] = (sumi[i - 1] + t) % mod;
sumG[i] = (sumG[i - 1] + quick\_pow(t, x) \* g[i] % mod) % mod;
}
}

void solve() {
cin >> T >> k >> x;
Mobi();
while(T--) {
int n;
cin >> n;
ll ans = 0;
for(int l = 1, r;l <= n; l = r + 1) {
r = min(n, n / (n / l));
ans = (ans + (sumG[r] - sumG[l - 1] + mod) % mod \* quick\_pow(sumi[n / l], x) % mod) % mod;
}

cout << ans << endl;
}
}

signed main() {
ios\_base::sync\_with\_stdio(false);
//cin.tie(nullptr);
//cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

huntian oy

$$求解\sum_{i=1}^n\sum_{j=1}^igcd(i^a-j^a,i^b-j^b)[gcd(i,j)=1]$$

$$gcd(i^a-j^a,i^b-j^b)=i^{gcd(a,b)}-j^{gcd(a,b)}$$

$题目说gcd(a,b)=1,所以那个复杂的式子直接变成i-j.$

$$\sum_{i=1}^n\sum_{j=1}^i(i-j)[gcd(i,j)=1]$$

$$\sum_{i=1}^n\sum_{j=1}^ii[gcd(i,j)=1]-\sum_{i=1}^n\sum_{j=1}^ij[gcd(i,j)=1]$$

$左边为与i互质的个数,右边为比i小并且互质的数字和,即为:$

$$\sum_{i=1}^ni\varphi (i)-\frac{\sum_{i=1}^ni\varphi (i)+1}{2}$$

$$\frac{\sum_{i=1}^ni\varphi (i)-1}{2}$$

$因为n有1e9,所以需要杜教筛优化。$

$$f*g(n)=\sum_{dn}f(d)g(\frac{n}{d})=\sum_{dn}id(d)*\varphi (d)*g(\frac{n}{d})$$

$很显然,设g=id,则:$

$$f*g(n)=\sum_{dn}\varphi (d)id(n)=n\sum_{dn}\varphi (d)=id(n)(\varphi *I)=id^2(n)=n^2$$

$$所以得:S(n)=\sum_{i=1}^n(f*g)i-\sum_{i=2}^ng(i)S(\frac{n}{i})=\frac{n(n+1)(2n+1)}{6}-\sum_{i=2}^niS(\frac{n}{i})$$

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
109
110
111
112
113
114
115
##include "bits/stdc++.h"
using namespace std;

typedef long long ll;
typedef pair<int, int> pii;
typedef pair<double, double> pdd;
typedef pair<ll, ll> pll;

##define endl "\n"
##define eb emplace\_back
##define mem(a, b) memset(a , b , sizeof(a))

const ll INF = 0x3f3f3f3f;
// const ll mod = 998244353;
const ll mod = 1e9 + 7;
const double eps = 1e-6;
const double PI = acos(-1);
const double R = 0.57721566490153286060651209;

template<typename T>
inline void read(T &a) {char c = getchar();T x = 0, f = 1;
while (!isdigit(c)) {if (c == '-')f = -1;c = getchar();}
while (isdigit(c)) {x = (x << 1) + (x << 3) + c - '0';c = getchar();}a = f \* x;
}

const int N = 1e6 + 10;

bool is\_prime[N];
int prime[N], cnt;
int phi[N];
ll sum\_i\_phi[N];

void init() {
phi[1] = 1;
for(int i = 2;i < N; i++) {
if(!is\_prime[i]) prime[++cnt] = i, phi[i] = i - 1;
for(int j = 1;j <= cnt && i \* prime[j] < N; j++) {
is\_prime[i \* prime[j]] = 1;
if(i % prime[j] == 0) {
phi[i \* prime[j]] = phi[i] \* prime[j];
break;
}
else phi[i \* prime[j]] = phi[i] \* phi[prime[j]];
}
}
for(int i = 1;i < N; i++) {
sum\_i\_phi[i] = (sum\_i\_phi[i - 1] + 1ll \* i \* phi[i] % mod) % mod;
}
}

ll quick\_pow(ll a, ll b) {
ll ans = 1;
while(b) {
if(b & 1) ans = ans \* a % mod;
a = a \* a % mod;
b >>= 1;
}
return ans % mod;
}

ll inv2;
ll inv6;

map<ll, ll> mp;

ll S(ll x) {
if(x < N) return sum\_i\_phi[x];
else if(mp[x]) return mp[x];
ll ans = x \* (x + 1) % mod \* (2 \* x + 1) % mod \* inv6 % mod;
for(int l = 2, r;l <= x; l = r + 1) {
r = min(x, x / (x / l));
ans = (ans - 1ll \* (r + l) % mod \* (r - l + 1) % mod \* inv2 % mod \* S(x / l)) % mod;
}
return mp[x] = ans;
}

void solve() {
init();
inv2 = quick\_pow(2, mod - 2);
inv6 = quick\_pow(6, mod - 2);
int \_; scanf("%d",&\_);
while(\_--) {
ll n, a, b; scanf("%lld%lld%lld",&n,&a,&b);
ll ans = (S(n) - 1) % mod;
ans = ans \* inv2 % mod;
printf("%lld\n",(ans % mod + mod) % mod);
}
}

signed main() {
ios\_base::sync\_with\_stdio(false);
// cin.tie(nullptr);
// cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

HDU 6428 Problem C. Calculate

$$求解\sum_{i=1}^A\sum_{j=1}^B\sum_{k=1}^C\phi (gcd(i,j^2,k^3))\;\;mod\;\;2^{30}$$

$首先要把\phi和gcd分开,即\phi(n)=\sum_{dn}(\phi*\mu)(d),证明如下:$ $$\phi(n)=(\mu*id)(n)$$

$$\sum_{dn}\mu(d)\frac{n}{d}$$

$$\sum_{dn}\mu(d)\sum_{d’\frac{n}{d}}\phi(d’)$$

$$\sum_{dd’n}\mu(d)\phi(d’)$$

$$\sum_{dn}\sum_{d’d}\mu(d’)\phi(\frac{d}{d’})$$

$$\sum_{dn}(\phi*\mu)(d)$$

$所以原式变为:$

$$\sum_{i=1}^A\sum_{j=1}^B\sum_{k=1}^C\sum_{di\;dj^2\;dk^3}(\phi*\mu)(d)$$

$$\sum_{d=1}^{A}(\phi*\mu)(d)\sum_{i=1\;di}^A\sum_{j=1\;dj^2}^B\sum_{k=1\;dk^3}^C1$$ \left \lceil \right \rceil $观察后面式子,对于一个x^k,若dx^k,先把d分解为\prod p_i^{a_i}.$ $则\prod p_i^{a_i}x^k,得\prod p_i^{\left \lceil \frac{a_i}{k} \right \rceil }x,所以设$

$$f_k(n)=\prod p_i^{\left \lceil \frac{a_i}{k} \right \rceil }$$

$则:$

$$ans=\sum_{d=1}^{A}(\phi * \mu)(d)\frac{A}{f_1(d)}\frac{B}{f_2(d)}\frac{C}{f_3(d)}$$

$对于f_k(n),类似分解n的形式。$ $分解质因子的过程中,记录质因子的指数,每次质因子+1时,若\%k=1时,说明该向上取整了,于是f_k(n)*=该质因子.$

$由于\phi和\mu都是积性函数,卷积之后还是积性函数,设g(d)=(\phi*\mu)(d),则$

$$\phi(n)=\sum_{dn}g(d)$$

$那么可以得:$

$$\phi(p^k)=\phi(p^{k-1})+g(p^k)$$

$$g(p^k)=\phi(p^k)-\phi(p^{k-1})$$

$$g(p^k)=(p-1)p^{k-1}-(p-1)p^{k-2}$$

$$g(p^k)=(p-1)^2p^{k-2}$$

$当我们欧拉筛的过程中:$

$$g(1)=1$$

$$g(p)=p-2$$

$$g(p^{k})=(p-1)^2p^{k-2}$$

$$g(p_1^{k_1}p_2^{k_2})=g(p_1^{k_1})g(p_2^{k_2})$$

$$……$$

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"
using namespace std;

typedef long long ll;

const ll mod = (1ll << 30);

const int N = 1e7 + 10;
bool is\_prime[N];
int prime[N], cnt;
int g[N];
int f1[N], f2[N], f3[N];
int deg[N];

void init() {
f1[1] = f2[1] = f3[1] = g[1] = 1;
for(int i = 2;i < N; i++) {
f1[i] = i;
if(!is\_prime[i]) prime[++cnt] = f2[i] = f3[i] = i, g[i] = i - 2, deg[i] = 1;
for(int j = 1;j <= cnt && i \* prime[j] < N; j++) {
int now = i \* prime[j];
is\_prime[now] = 1;
if(i % prime[j] == 0) {
deg[now] = deg[i] + 1;
int num = 1, tmp = i;
while(num <= 3 && tmp % prime[j] == 0) num++, tmp /= prime[j];
if(num == 1) g[now] = g[i] \* g[prime[j]];
else if(num == 2) g[now] = g[i / prime[j]] \* (prime[j] - 1) \* (prime[j] - 1);
else g[now] = g[i] \* prime[j];
f2[now] = f2[i] \* (deg[now] % 2 == 1 ? prime[j] : 1);
f3[now] = f3[i] \* (deg[now] % 3 == 1 ? prime[j] : 1);
break;
}
else {
deg[now] = 1;
f2[now] = f2[i] \* f2[prime[j]];
f3[now] = f3[i] \* f3[prime[j]];
g[now] = g[i] \* g[prime[j]];
}
}
}
}

void solve() {
init();
int \_; scanf("%d",&\_);
while(\_--) {
int A, B, C; scanf("%d%d%d",&A,&B,&C);
ll ans = 0;
for(int d = 1;d <= A; d++) {
ans = (ans + 1ll \* g[d] \* (A / f1[d]) % mod \* (B / f2[d]) % mod \* (C / f3[d]) % mod) % mod;
}
printf("%lld\n",(ans % mod + mod) % mod);
}
}

signed main() {
ios\_base::sync\_with\_stdio(false);
// cin.tie(nullptr);
// cout.tie(nullptr);
##ifdef FZT\_ACM\_LOCAL
freopen("in.txt", "r", stdin);
freopen("out.txt", "w", stdout);
signed test\_index\_for\_debug = 1;
char acm\_local\_for\_debug = 0;
do {
if (acm\_local\_for\_debug == '$') exit(0);
if (test\_index\_for\_debug > 20)
throw runtime\_error("Check the stdin!!!");
auto start\_clock\_for\_debug = clock();
solve();
auto end\_clock\_for\_debug = clock();
cout << "Test " << test\_index\_for\_debug << " successful" << endl;
cerr << "Test " << test\_index\_for\_debug++ << " Run Time: "
<< double(end\_clock\_for\_debug - start\_clock\_for\_debug) / CLOCKS\_PER\_SEC << "s" << endl;
cout << "--------------------------------------------------" << endl;
} while (cin >> acm\_local\_for\_debug && cin.putback(acm\_local\_for\_debug));
##else
solve();
##endif
return 0;
}

本文作者:jujimeizuo
本文地址https://blog.jujimeizuo.cn/2020/08/10/mobius-inversion/
本博客所有文章除特别声明外,均采用 CC BY-SA 3.0 协议。转载请注明出处!