题解 P5176 【公约数】

Nemlit

2019-09-27 22:28:20

Solution

~~看到这题目就不想做了系列,出题人是不是都不知道杜教筛是什么东西啊,他家杜教筛可以预处理优化到~~$O(n^\frac{2}{3})$ ~~先吓唬你一下,我们要求:~~ $$\sum_{i=1}^n\sum_{j=1}^m\sum_{k=1}^p\gcd(i\cdot j,i\cdot k,j\cdot k)\times \gcd(i,j,k)\times \left(\frac{\gcd(i,j)}{\gcd(i,k)\times \gcd(j,k)}+\frac{\gcd(i,k)}{\gcd(i,j)\times \gcd(j,k)}+\frac{\gcd(j,k)}{\gcd(i,j)\times \gcd(i,k)}\right)$$ 我们令$x = gcd(i, j), y = gcd(i, k), z = gcd(j, k)$ 于是我们有: $$\sum_{i=1}^n\sum_{j=1}^m\sum_{k=1}^p\gcd(i\cdot j,i\cdot k,j\cdot k)\times \gcd(i,j,k)\times \left(\frac{x}{y\times z}+\frac{y}{x\times z}+\frac{z}{x\times y}\right)$$ $$\sum_{i=1}^n\sum_{j=1}^m\sum_{k=1}^p\gcd(i\cdot j,i\cdot k,j\cdot k)\times \gcd(i,j,k) \times \frac{x ^ 2 + y ^ 2 + z ^ 2}{x \times y \times z}$$ 展开第一陀恶心人的柿子,可以得到: $$gcd(i\cdot j,i\cdot k,j\cdot k) = \frac{x \times y\times z}{gcd(i, j, k)}$$ (证明可以通过**唯一分解定理**展开) 于是我们惊喜地发现,我们只需要求$($还原$x, y, z)$: $$\sum_{i=1}^n\sum_{j=1}^m\sum_{k=1}^p gcd(i, j)^2 +gcd(j, k) ^ 2 + gcd(i, k) ^ 2$$ 我们把柿子分开写,得到: $$\sum_{i=1}^n\sum_{j=1}^m\sum_{k=1}^p gcd(i, j)^2+\sum_{i=1}^n\sum_{j=1}^m\sum_{k=1}^pgcd(j, k) ^ 2 +\sum_{i=1}^n\sum_{j=1}^m\sum_{k=1}^p gcd(i, k) ^ 2$$ 因为$gcd(i, j)$与k无关(其他同理),于是我们可以提出k,变成p个$\sum_{i=1}^n\sum_{j=1}^m gcd(i, j)^2$相加,故有: $$p \times \sum_{i=1}^n\sum_{j=1}^m gcd(i, j)^2+n \times \sum_{j=1}^m\sum_{k=1}^pgcd(j, k) ^ 2 +m \times \sum_{i=1}^n\sum_{k=1}^p gcd(i, k) ^ 2$$ 设$F(n, m) = \sum_{i=1}^n\sum_{j=1}^m gcd(i, j)^2$ 原式$=p \times F(n, m) + n \times F(p, m) + m \times F(n, p)$ 所以问题转化为如何快速的求出$F(n, m)$(不妨设$n ≤ m$) $$F(n, m) = \sum_{d = 1}^n\sum_{i=1}^n\sum_{j=1}^m d^2 \times [gcd(i, j) == d]$$ $i, j$同时除以d得: $$F(n, m) = \sum_{d = 1} ^ n\sum^{n / d}_{i = 1}\sum^{m / d}_{j = 1} d ^ 2 \times [gcd(i, j) == 1]$$ $$F(n, m) = \sum_{d = 1} ^ nd ^ 2\sum^{n / d}_{i = 1}\sum^{m / d}_{j = 1}\sum_{k|gcd(i, j)} \mu(k)$$ 提前枚举K得: $$F(n, m) = \sum_{d = 1} ^ n d ^ 2\sum_{k = 1}^{n / d}\mu(k)\sum^{n / (d * k)}_{i = 1}\sum^{m / (d * k)}_{j = 1} $$ $$F(n, m) = \sum_{d = 1} ^ n d ^ 2\sum_{k = 1}^{n / d}\mu(k)\left \lfloor \frac{n}{d * k} \right \rfloor \times \left \lfloor \frac{m}{d * k} \right \rfloor$$ 跟据反演套路,我们令$T = d * k$ $$F(n, m) = \sum_{T = 1}^n \left \lfloor \frac{n}{T} \right \rfloor \times \left \lfloor \frac{m}{T} \right \rfloor\sum_{d | T} d ^ 2\times \mu(\frac{T}{d})$$ 令$g(x) = x * x, mu(x) = \mu(x)$ 于是$F(n, m) = \sum_{T = 1}^n \left \lfloor \frac{n}{T} \right \rfloor \times \left \lfloor \frac{m}{T} \right \rfloor\times(g(x) \times mu(x))$ 发现后面那一堆就是一个卷积的形式,而且是两个积性函数的卷积,于是我们可以用线性筛筛出,筛法见[于神之怒加强版](https://www.luogu.org/problem/P4449) 总结一下求法: 我们记$prim[i]$为第i个质数,$low[i]$为i质因数分解后最小的质因子的指数次幂,$f[x]$为$g(x)$卷$mu(x)$,$mu[x]$为$\mu(x)$ 类似于线性筛,我们分三种情况讨论: $case1:$ $i$是质数 我们可以直接算出来:$low[i] = i,\;\ f[i] = i \times i - 1,\;\ mu[i] = -1$ $case2:i \% prim[j] != 0$:这种情况是满足积性函数定义的,即$i, prim[j]$互质,令$p = i \times prim[j]$,所以我们有$low[p] = prim[j],\;\ mu[p] = -mu[i],\;\ f[i] = f[i] \times f[prim[j]]$ $case3:i \% prim[j] == 0$:这种情况不满足积性函数的定义,于是我们考虑用我们刚刚维护的$low[i]$ 考虑到$i / low[i],prim[j]$是互质的,那我们可不可以类比$case2$通过$f[i / low[i]], f[prim[j]]$来推出$f[p]$呢? 假设$a, b$互质,我们有:$f[a \times b] = f[a] \times f[b]$ 所以我们也有:$f[p] = f[i] \times f[prim[j]] = f[i\ /\ low[i]] \times f[prim[j] \times low[i]]$ 根据定义,$i\ /\ low[i]$没有$prim[j]$这个约数,所以$i / low[i],\;\ prim[j] \times low[i]$互质 i肯定有$prim[j]$这个质因数,而且根据线性筛的性质,$prim[j]$肯定是i的最小质因数,故我们有$mu[p] = 0,\;\ low[p] = low[i] \times prim[j]$ 如果$low[i] < i$,证明$low[i] \times prim[j] <= i$,故我们可以直接用$f[p] = f[i\ /\ low[i]] \times f[prim[j] \times low[i]]$ 如果$low[i] == i$,证明i以前肯定没有被筛过,若$low[i] = prim[j] ^ k$,则k肯定为质数,那么我们就可以直接用:$f[p] = (i * prim[j]) ^ 2 - i ^ 2$ ### $ps:$ 写完后看题解才发现:这道题根于神之怒不一样,并不需要记录$low[i]$,对于$case3$,$f[p] = f[i] \times prim[j] ^2$即可 # $Code:$ ``` #include<bits/stdc++.h> using namespace std; #define il inline #define re register il int read() { re int x = 0, f = 1; re char c = getchar(); while(c < '0' || c > '9') { if(c == '-') f = -1; c = getchar();} while(c >= '0' && c <= '9') x = x * 10 + c - 48, c = getchar(); return x * f; } #define rep(i, s, t) for(re int i = s; i <= t; ++ i) #define mod 1000000007 #define maxn 20000007 il int Pow(int x) {return 1ll * x * x % mod;} int cnt, prim[maxn], f[maxn], mu[maxn], low[maxn]; bool vis[maxn]; il void init() { f[1] = mu[1] = 1; rep(i, 2, 2e7) { if(!vis[i]) f[i] = Pow(i) - 1, mu[i] = -1, low[i] = i, prim[++ cnt] = i; for(re int j = 1; j <= cnt; ++ j) { if(i * prim[j] > 2e7) break; int p = i * prim[j]; vis[p] = 1; if(i % prim[j] == 0) { low[p] = low[i] * prim[j]; if(low[i] == i) f[p] = (Pow(1ll * i * prim[j] % mod) - Pow(i) + mod) % mod; else f[p] = 1ll * f[i / low[i]] * f[low[i] * prim[j]] % mod; break; } mu[p] = -mu[i], low[p] = prim[j], f[p] = 1ll * f[i] * f[prim[j]] % mod; } } rep(i, 1, 2e7) f[i] = (f[i - 1] + f[i]) % mod; } il int F(int n, int m) { if(n > m) swap(n, m); int ans = 0; for(re int l = 1, r; l <= n; l = r + 1) { r = min(n / (n / l), m / (m / l)); ans = (ans + 1ll * (1ll * (n / l) * (m / l) % mod) * (f[r] - f[l - 1]) % mod) % mod; } return (ans + mod) % mod; } il int get(int n, int m, int p) { return ((1ll * p * F(n, m) % mod + 1ll * n * F(p, m) % mod) % mod + 1ll * m * F(n, p)) % mod; } il void outit() { for(re int i = 1, T = read(); i <= T; ++ i) printf("%d\n", get(read(), read(), read())); } int main() { init(), outit(); return 0; } ```