P4619 [SDOI2018]旧试题
题目描述
求下面表达式的值
[largesum_{i=1}^{A}sum_{j=1}^{B}sum_{k=1}^{C}sigma_{0}(ijk) mod(10^9+7)
]
数据范围
(1 le A,B,C le 100000)
解题思路
这题非常困难是因为推着推着式子变成了图论问题？
让我们先来推一波式子把
要用到的结论
[large sigma_{0}(ijk)=sum_{x|i}sum_{y|j}sum_{z|k}[gcd(x,y)=1][gcd(y,z)=1][gcd(x,z)=1]\
[x=1] = sum_{d|x}mu(d)
]
第一个式子的理解
单独考虑一个质因子 P，i 中最大有 (P^a)，j 中最大有 (P^b)，k 中最大有 (P^c)，那么 (ijk) 中有 (P^{a+b+c})， 现在我们尝试建立起映射，选因数时，如果 P 选了小于等于 a 个，那么全在 i 中选，显然满足全部互质的条件，大于 a 个但小于等于 a + b 个，那么优先在 i 中选，剩下的在 j 中选，然后将 i 清零，大于 a + b 也同理，也就是选质数时尽量靠前选，但我们只让最后一个被选到的有，前面的设为 0，这样可以保证互质，且是一种双射关系
第二个式子的理解
裸的莫比乌斯反演
正式开始推式子
[largesum_{i=1}^{A}sum_{j=1}^{B}sum_{k=1}^{C}sigma_{0}(ijk)\
large=sum_{i=1}^{A}sum_{j=1}^{B}sum_{k=1}^{C}sum_{x|i}sum_{y|j}sum_{z|k}[gcd(x,y)=1][gcd(y,z)=1][gcd(x,z)=1]\
large=sum_{x=1}^{A}sum_{y=1}^{B}sum_{z=1}^{C}[gcd(x,y)=1][gcd(y,z)=1][gcd(x,z)=1]sum_{i=1}^{frac Ax}sum_{j=1}^{frac By}sum_{k=1}^{frac Cz}\
large=sum_{x=1}^{A}sum_{y=1}^{B}sum_{z=1}^{C}[gcd(x,y)=1][gcd(y,z)=1][gcd(x,z)=1]lfloorfrac Ax
floorlfloorfrac By
floorlfloorfrac Cz
floor\
large=sum_{x=1}^{A}sum_{y=1}^{B}sum_{z=1}^{C}sum_{u|x,u|y}mu(u)sum_{v|y,v|z}mu(v)sum_{w|x,w|z}mu(w)lfloorfrac Ax
floorlfloorfrac By
floorlfloorfrac Cz
floor\
large=sum_{u=1}^{min(A,B)}sum_{v=1}^{min(B,C)}sum_{w=1}^{min(A,C)}mu(u)mu(v)mu(w)sum_{x=lcm(u,w)}lfloorfrac Ax
floorsum_{y=lcm(u,v)}lfloorfrac By
floorsum_{z=lcm(v,w)} lfloorfrac Cz
floor\
large=sum_{u=1}^{min(A,B)}sum_{v=1}^{min(B,C)}sum_{w=1}^{min(A,C)}mu(u)mu(v)mu(w)F(lcm(u,w),A)F(lcm(u,v),B)F(lcm(v,w),C)
]
然后就是超级神奇的操作
我们把枚举的三元组抽象为一个三元环 (x, y, z)，两个点间连边权为 lcm 的边，这个三元环对答案有贡献当且仅当
- (mu(x) eq0 ~~mu(y) eq0~mu(z) eq0)
- (lcm(x, y, z) le max(A, B, C))
你觉得它可能会快，但你认为它还是暴力，事实上它跑的飞快（比纯暴力来说
枚举三元环肯定不是暴力枚举，容易发现 lcm 的 (mu) 也不为零，我们可以直接枚举 lcm，将其质因数分解，根据质因数枚举其中一个端点，然后暴力枚举另一个端点，这样枚举的每一个都是合法的点对，时间复杂度 (Theta(m))
打表发现有大概 70 多万条边，暴力枚举肯定不行，用上黑科技 三元环计数 可以做到 (Theta(msqrt m))，卡卡常勉强跑过。
一发过了，但除法真的比乘法慢一大堆，请谨慎使用除法！！
#pragma GCC optimize(3, "inline")
#include <queue>
#include <vector>
#include <iostream>
#include <cstdio>
#include <cstring>
#include <algorithm>
#define MP make_pair
#define ll long long
#define fi first
#define se second
using namespace std;
template <typename T>
void read(T &x) {
x = 0; bool f = 0;
char c = getchar();
for (;!isdigit(c);c=getchar()) if (c=='-') f=1;
for (;isdigit(c);c=getchar()) x=x*10+(c^48);
if (f) x=-x;
}
template<typename F>
inline void write(F x)
{
static short st[30];short tp=0;
if(x<0) putchar('-'),x=-x;
do st[++tp]=x%10,x/=10; while(x);
while(tp) putchar('0'|st[tp--]);
putchar('
');
}
template <typename T>
inline void Mx(T &x, T y) { x < y && (x = y); }
template <typename T>
inline void Mn(T &x, T y) { x > y && (x = y); }
#define re register
const int Maxn = 300000, N = Maxn + 5;
const int P = 1e9+7;
ll f[N], mu[N], e[N], prime[N], p[N], tot;
void prework(void) {
mu[1] = 1; f[1] = 1;
for (re int i = 2;i <= Maxn; i++) {
if (!e[i]) prime[++tot] = e[i] = i, mu[i] = -1, f[i] = p[i] = 2;
for (re int j = 1; prime[j] * i <= Maxn && j <= tot; j++) {
int t = prime[j] * i; e[t] = prime[j];
if (prime[j] == e[i]) { p[t] = p[i] + 1, f[t] = f[i] / p[i] * p[t]; break; }
f[t] = f[i] * (p[t] = 2), mu[t] = mu[i] * -1;
}
// cout << f[i] << endl;
f[i] = f[i] + f[i-1]; if (f[i] >= P) f[i] -= P;
}
}
struct node {
int x, y, lcm;
}ed[N*20];
ll ans;
ll vis[N], lm[N], a[N], b[N], c[N], A, B, C;
ll st[N], deg[N], lim, tp, T;
vector<pair<ll, ll> > v[N];
inline void clear(void) {
lim = max(max(A, B), C), tot = ans = 0;
for (int i = 1;i <= lim; i++)
a[i] = f[A/i], b[i] = f[B/i], c[i] = f[C/i], v[i].clear();
}
inline ll calc(int x, int y, int z) {
return a[x] * b[y] * c[z];
}
inline void build(void) {
// cout << calc(2, 1, 2) + calc(2, 2, 1) + calc(1, 2, 2) << endl;
memset(deg, 0, lim * 4 + 200);
for (int i = 1;i <= lim; i++) {
if (!mu[i]) continue; int x = i; tp = 0;
// if (i % 1000 == 0) cout << i << endl;
while (x != 1) st[++tp] = e[x], x /= e[x];
for (re int j = 0; j < (1 << tp); j++) {
int u = 1, t, v;
for (re int k = 1;k <= tp; k++)
if (j & (1 << (k - 1))) u *= st[k];
t = i / u;
for (int s = j; ; s = (s - 1) & j) {
v = t;
for (int p = 0;p < tp; p++)
if (s >> p & 1) v *= st[p + 1];
// cout << u << ' ' << v << endl;
if (u > v) {
ed[++tot] = (node) {u, v, i}; deg[u]++, deg[v]++;
// cout << u << ' ' << i << endl;
// cout << calc(u, i, i) + calc(i, u, i) + calc(i, i, u) << endl;
ans += mu[v] * (calc(u, i, i) + calc(i, u, i) + calc(i, i, u));
ans += mu[u] * (calc(v, i, i) + calc(i, v, i) + calc(i, i, v));
}
if (!s) break;
}
}
ans += mu[i] * calc(i, i, i);
}
for (int i = 1;i <= tot; i++) {
int x = ed[i].x, y = ed[i].y;
if (deg[x] > deg[y] || (deg[x] == deg[y] && x > y)) swap(x, y);
v[x].emplace_back(y, ed[i].lcm);
}
}
void work(void) {
memset(vis, 0, lim * 4 + 200);
// cout << calc(10, 5, 2) + calc(10, 2, 5) + calc(2, 5, 10) + calc(2, 10, 5) << endl;
for (int i = 1;i <= lim; i++) {
for (pair<int, int> to: v[i])
vis[to.fi] = i, lm[to.fi] = to.se;
for (pair<int, int> to: v[i]) {
for (pair<int, int> To: v[to.fi]) {
if (vis[To.fi] != i) continue;
// int pre = ans;
// cout << lm[To.fi] << ' ' << to.se << ' ' << To.se << endl;
ans += mu[i] * mu[to.fi] * mu[To.fi] *
( calc(lm[To.fi], to.se, To.se)
+ calc(lm[To.fi], To.se, to.se)
+ calc(to.se, lm[To.fi], To.se)
+ calc(to.se, To.se, lm[To.fi])
+ calc(To.se, lm[To.fi], to.se)
+ calc(To.se, to.se, lm[To.fi]) );
// cout << ans - pre << endl;
}
}
ans %= P;
}
ans = (ans % P + P) % P;
}
int main() {
for (prework(), read(T); T; T--) {
read(A), read(B), read(C);
clear(), build(), work(), write(ans);
}
return 0;
}
/*
1
100000 100000 100000
*/