【题解】旧试题

来自神仙 11Dimensions 的神仙做法……

题目

原题链接

求:
$$ \sum\limits_{i=1}^A \sum\limits_{j=1}^B \sum\limits_{k=1}^C d(ijk) $$
其中 $A,B,C ⩽ 200005$。

分析

记$x⊥y$表示$(x,y)=1$。
结论:

$$d(ijk)=\sum\limits_{x|i}\sum\limits_{y|i}\sum\limits_{z|i}[x⊥y][x⊥z][y⊥z]$$

其中[]表示艾弗森括号,当且仅当括号内命题为真时取值为1,否则为0。
在本篇题解中,这个括号的意义等价于单位函数,即$\epsilon((a,b))=[(a,b)=1]=[a⊥b]$。
那么:

$$\begin{aligned}& \sum\limits_{i=1}^A \sum\limits_{j=1}^B \sum\limits_{k=1}^C d(xyz) \newline
=& \sum\limits_{i=1}^A \sum\limits_{j=1}^B\sum\limits_{k=1}^C\sum\limits_{x|i}\sum\limits_{y|i}\sum\limits_{z|i}[x⊥y][x⊥z][y⊥z] \newline
=& \sum\limits_{x=1}^A \sum\limits_{y=1}^B \sum\limits_{z=1}^C[x⊥y][x⊥z][y⊥z] \lfloor\frac{A}{x}\rfloor \lfloor\frac{B}{y}\rfloor \lfloor\frac{C}{z}\rfloor \end{aligned}$$

从三个单位函数里面任选一个反演,利用$\epsilon=\mu*1$:

$$\begin{aligned}& \sum\limits_{x=1}^A \sum\limits_{y=1}^B \sum\limits_{z=1}^C[x⊥y][x⊥z][y⊥z] \lfloor\frac{A}{x}\rfloor \lfloor\frac{B}{y}\rfloor \lfloor\frac{C}{z}\rfloor \newline
=& \sum\limits_{x=1}^A \sum\limits_{y=1}^B \sum\limits_{z=1}^C (\sum\limits_{d|x,d|y}\mu(d)) [x⊥z][y⊥z] \lfloor\frac{A}{x}\rfloor \lfloor\frac{B}{y}\rfloor \lfloor\frac{C}{z}\rfloor \newline
=& \sum\limits_{z=1}^C\lfloor\frac{C}{z}\rfloor\sum\limits_{d=1}^{\min(A,B)}\mu(d)\sum\limits_{k_1=1}^{\lfloor\frac{A}{d}\rfloor}\sum\limits_{k_2=1}^{\lfloor\frac{B}{d}\rfloor}[k_1d⊥z][k_2d⊥z]
\lfloor\frac{A}{k_1d}\rfloor \lfloor\frac{B}{k_2d}\rfloor\end{aligned}$$

依据$[ab⊥c]\iff[a⊥c][b⊥c]$,整理得到:

$$\sum\limits_{z=1}^C\lfloor\frac{C}{z}\rfloor\sum\limits_{d=1}^{\min(A,B)}\mu(d)[d⊥z]
(\sum\limits_{x=1}^{\lfloor\frac{A}{d}\rfloor}[x⊥z]\lfloor\frac{A}{xd}\rfloor)
(\sum\limits_{y=1}^{\lfloor\frac{B}{d}\rfloor}
[y⊥z] \lfloor\frac{B}{yd}\rfloor)$$

记:

$$\begin{aligned} g(n,x) = \sum\limits_{i=1}^n\mu(i)[i⊥x] \newline
f(n,x) = \sum\limits_{i=1}^n\lfloor\frac{n}{i}\rfloor[i⊥x] \end{aligned}$$

答案变成:

$$\begin{aligned}& \sum\limits_{z=1}^C\lfloor\frac{C}{z}\rfloor\sum\limits_{d=1}^{\min(A,B)}\mu(d)[d⊥z]f(\lfloor\frac{A}{d}\rfloor, z)f(\lfloor\frac{B}{d}\rfloor, z) \newline
=& \sum\lfloor\frac{C}{z}\rfloor\sum(g(r,z) - g(l - 1, z))f(\lfloor\frac{A}{d}\rfloor, z)f(\lfloor\frac{B}{d}\rfloor, z)\end{aligned}$$

即对$\mu()$做前缀和然后对后面的$f()$分段,其中$[l,r]$表示整除分段的一段区间。

假设$O(n)$枚举$z$,那么求出后面的$\sum$的值是$O(\sqrt{n})$的。但是发现$z$是变化的,当$z$变化时暴力维护$f(),g()$是$O(n^2)$的,这成为了代码复杂度的瓶颈。如何解决?
能不能减少更新$f(),g()$的次数?不难发现$z$在函数中发挥作用的地方是判断一个数与其互质。考虑将$z$质因数分解,容易看出一个数与其互质仅和$z$的质因子的种类有关,而与质因子的幂次无关。
记:

$$lw(z) = \prod\limits_{i=1}^np_i, \text{where }z = \prod\limits_{i=1}^np_i^{\alpha_i}$$

那么我们只需要考虑$z\in [1,C]$的所有$lw(z)$值即可(即所有无平方因子的数),它们共用一套$f(),g()$的函数值。通过dfs暴力生成无平方因子数,我们可以把它们一并更新。

但是这远远不够,考虑通过递推来维护$f(),g()$。这个套路参考NOI2016 循环之美
考虑在$z$中删去其一个质因子$x$,即$f(n,z)\to f(n,z/x)$,将出现何种变化?应当有一部分多加了,要减去:

$$\begin{aligned}f(n,z) &= \sum\limits_{i=1}^n\lfloor\frac{n}{i}\rfloor[i⊥z/x] -
\sum\limits_{i=1}^n\lfloor\frac{n}{i}\rfloor[i⊥z/x][x|i] \newline
&= f(n, z / x) - \sum\limits_{k=1}^{\lfloor\frac{n}{x}\rfloor}[kx⊥z/x]\lfloor\frac{n}{kx}\rfloor \newline
&= f(n, z / x) - [x⊥z/x]\sum\limits_{k=1}^{\lfloor\frac{n}{x}\rfloor}[k⊥z/x]\lfloor\frac{n}{kx}\rfloor\newline
&= f(n, z / x) - f(\lfloor\frac{n}{x}\rfloor, z/x)\end{aligned}$$

对$g()$的推导也同理,利用$\mu(ab)=\mu(a)\mu(b)[a⊥b]$,可以得出:

$$\begin{aligned}g(n,z) &= \sum\limits_{i=1}^n \mu(i) [i⊥z/x] -
\sum\limits_{i=1}^n \mu(i) [i⊥z/x][x|i] \newline
&= g(n, z / x) - \sum\limits_{k=1}^{\lfloor\frac{n}{x}\rfloor}[kx⊥z/x] \mu(kx) \newline
&= g(n, z / x) - \mu(x)\sum\limits_{k=1}^{\lfloor\frac{n}{x}\rfloor}\mu(k)[k⊥x][k⊥z/x] \newline
&= g(n, z / x) + \sum\limits_{k=1}^{\lfloor\frac{n}{x}\rfloor}\mu(k)[k⊥x(z/x)] \newline
&= g(n, z / x) + g(\lfloor\frac{n}{x}\rfloor, z)\end{aligned}$$

即:

$$\begin{aligned} g(n,z) = g(n, z / x) + g(\lfloor\frac{n}{x}\rfloor, z) \newline
f(n,z) = f(n, z / x) - f(\lfloor\frac{n}{x}\rfloor, z/x) \end{aligned}$$

那么就可以递推了,但是直接更新是$O(n)$的。容易发现整除分段并不会用到所用的$f(),g()$值,所以我们边做分段边更新就好了,这是$O(\sqrt{n})$的。
空间复杂度呢?f[2e5][2e5],g[2e5][2e5]看似存不下来,但是容易发现后面一维是不连续使用的,那么将其离散化,考虑在dfs构造无平方因子数的时候,下层状态的转移依赖于上层状态,而dfs树的深度是$O(\log_2n)$的,所以开f[2e5][10]即可。

总复杂度?粗略估计,爆搜的复杂度是$O(2^{\log_2n})$即$O(n)$的,后面的求和通过整除分段可以$O(\sqrt{n})$得出,那么总复杂度是$O(n\sqrt{n})$的。但实际上复杂度要小很多,也就是说如果常数写得好那么它可以跑的飞快(预处理整除分段的端点、cache-friendly之类的),但是本人代码常数没那么小,最慢一个点大约4s?
比三元环是好得多了。

Orz 11Dimensions

代码

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<iostream>
#include<cstdio>
#include<cstring>
#include<string>
using namespace std;

#define LL long long

const int CN = 2e5+5;
const LL P = 1e9+7;

LL read(){
LL s = 0,ne = 1; char c = getchar();
for(;c < '0' || c > '9'; c = getchar()) if(c == '-') ne = -1;
for(;c >= '0' && c <= '9';c = getchar()) s = (s << 1) + (s << 3) + c - '0';
return s * ne;
}

int p[CN],mu[CN],lw[CN]; LL d[CN]; bool np[CN];
void sieve(int n){
np[1] = true; mu[1] = 1, d[1] = 1, lw[1] = 1;
for(int i = 2;i <= n; i++){
if(!np[i]) p[ ++p[0] ] = i, d[i] = 2, lw[i] = i, mu[i] = -1;
for(int j = 1;j <= p[0] && i * p[j] <= n; j++){
int x = i * p[j]; np[x] = true;
if(i % p[j]) d[x] = d[i] << 1, lw[x] = lw[i] * p[j], mu[x] = -mu[i];
else {d[x] = (d[i] << 1) - d[i / p[j]], lw[x] = lw[i]; break;}
}
}
for(int i = 1;i <= n;i++) d[i] += d[i - 1],mu[i] += mu[i - 1];
}

//
int t,A,B,C;
LL f[CN][10],g[CN][10],s[CN],ans;
void init(){
g[0][1] = f[0][1] = 0;
for(int l = 1;l <= A; l++){
int r = min(A / (A / l), B / (B / l));
g[r][1] = mu[r]; f[A / l][1] = d[A / l]; f[B / l][1] = d[B / l];
l = r;
}

memset(s, 0 , sizeof(s));
for(int z = 1;z <= C; z++) s[ lw[z] ] += 1ll * C / z;

ans = 0;
}
void upd(int x,int k){
for(int l = 1;l <= A; l++){
int r = min(A / (A / l), B / (B / l));
g[r][k] = g[r][k - 1] + g[r / x][k];
f[A / l][k] = f[A / l][k - 1] - f[(A / l) / x][k - 1];
f[B / l][k] = f[B / l][k - 1] - f[(B / l) / x][k - 1];
l = r;
}
}
void dfs(int z0, int u, int k){
LL cur = 0;
for(int l = 1;l <= A; l++){
int r = min(A / (A / l), B / (B / l));
cur += (g[r][k] - g[l - 1][k]) * f[A / l][k] * f[B / l][k];
cur = (cur % P + P) % P;
l = r;
}
ans += (s[z0] * cur) % P; ans %= P;
for(int v = u; 1ll * p[v] * z0 <= C; v++) upd(p[v], k + 1), dfs(z0 * p[v], v + 1, k + 1);
}

int main()
{
freopen("_in.in", "r", stdin);

t = read(); sieve(2e5);
while(t--){
A = read(), B = read(), C = read();
if(A > B) swap(A, B);

init(); dfs(1, 1, 1);
printf("%lld\n", ans);
}
}

$$ - - - - \mathcal{End} - - - - $$


评论

Your browser is out-of-date!

Update your browser to view this website correctly. Update my browser now

×