简单反演+神仙优化题。
首先我们知道
( 是除数函数,相当于题目描述中的 )
只要将 素因子分解后观察即可得到上面的结论。它的好处是出现了我们喜欢的
就可以莫比乌斯反演了
设 ,可以 预处理出 所有值
然而痛苦的地方就在于这样推并不能降复杂度,目前为止我们的复杂度依然是 。事实上可以把它推到 的,那样做可以过cf原题但是并不利于我们这题的后续优化。
但是也许我们可以不 地枚举 三元组,因为事实上这样枚举的过程中有很多情况的贡献是 :
-
如果 中有一个是 ,那么贡献是 。这样的情况还挺多的。
-
如果 ,那么 。也就是说 或者 或者 都会导致贡献是 。
所以我们有一种做法就是枚举最大公约数 ,枚举 并注意满足条件 ,连边 ,建出一张图。这张图的三元环刚好对应到 三元组,三元环计数的同时统计答案即可。
至于这样做的复杂度为什么是对的,我们已经发现有相当多的点对 并不满足上述条件(意味着它们产生的贡献是 ),除去这些之后剩下的点对连出的边不会太多。正如shadowice大佬已经说了他实测下来这个边数只有七十多万。那么只要你的三元环计数速度合格就好。
关于三元环计数,可以参考iki9的这篇博客,我们可以在 的时间内完成三元环计数。
但是注意到上述连边没有考虑到自环。也就是说如果三元组 中存在两个数相等或者三个数都相等,我们并没有将答案计算在内。因此我们需要单独算。这个很方便:三个数都相等的情况直接枚举这个数是几,把答案加起来就好;恰好有两个数相等的情况就在我们枚举 连边的过程中解决,只要顺手把 和 的答案都加起来即可。
在做的过程中时刻注意剪枝,只要出现 值为 就跳过。另外可以顺手加一些循环展开、register之类的卡卡常。还有大家应该都注意到的cache miss问题,我们选择使用vector而不是邻接表存图以访问连续的内存来加速。
还有一个小trick:我们统计的是 的约数个数之和,一个数的约数个数不会很多, 又都小于等于 ,估计最终的答案不会太大。事实上答案是不会超过 的,因此我们在统计的过程中完全不需要取模,只要在最后输出的时候模一下即可。
给一份稍微可读一点的代码
#include <cctype>
#include <cstdio>
#include <climits>
#include <algorithm>
#include <vector>
template <typename T> inline void read(T& x) {
int f = 0, c = getchar(); x = 0;
while (!isdigit(c)) f |= c == '-', c = getchar();
while (isdigit(c)) x = x * 10 + c - 48, c = getchar();
if (f) x = -x;
}
template <typename T> void write(T x) {
if (x < 0) putchar('-'), x = -x;
if (x > 9) write(x / 10);
putchar(x % 10 + 48);
}
template <typename T> inline void writeln(T x) { write(x); puts(""); }
template <typename T> inline bool chkmin(T &x, const T &y) { return y < x ? (x = y, 1) : 0; }
template <typename T> inline bool chkmax(T &x, const T &y) { return x > y ? (x = y, 1) : 0; }
#if __cplusplus >= 201103L
#define reg
#else
#define reg register
#endif
typedef long long LL;
const int maxn = 1e5 + 207;
const LL mod = 1e9 + 7;
int mu[maxn], pri[maxn];
bool ff[maxn];
LL fa[maxn], fb[maxn], fc[maxn];
int a, b, c, mx, mn;
struct Edge {
int u, v, w;
Edge(int x, int y, int z = 0) : u(x), v(y), w(z) {}
Edge() : u(0), v(0), w(0) {}
};
struct Node {
int to, w;
Node(int x, int y = 0) : to(x), w(y) {}
Node() : to(0), w(0) {}
};
std::vector<Node> G[maxn];
int deg[maxn];
Edge es[maxn << 4];
int vis[maxn];
int etot;
int gcd(int a, int b) { return b ? gcd(b, a % b) : a; }
inline void sieve(int n) {
mu[1] = 1;
for (reg int i = 2; i <= n; ++i) {
if (!ff[i]) { pri[++pri[0]] = i; mu[i] = -1; }
for (reg int j = 1; j <= pri[0]; ++j) {
int x = i * pri[j];
if (x > n) break;
ff[x] = 1;
if (i % pri[j]) mu[x] = -mu[i];
else { mu[x] = 0; break; }
}
}
}
int main() {
int T; read(T);
sieve(1e5);
while (T--) {
read(a); read(b); read(c);
mx = std::max(std::max(a, b), c);
mn = std::min(std::min(a, b), c);
// 预处理 f 函数值
for (reg int i = 1; i <= mx; i += 4) {
fa[i] = 0; fa[i + 1] = 0; fa[i + 2] = 0; fa[i + 3] = 0;
fb[i] = 0; fb[i + 1] = 0; fb[i + 2] = 0; fb[i + 3] = 0;
fc[i] = 0; fc[i + 1] = 0; fc[i + 2] = 0; fc[i + 3] = 0;
}
for (reg int i = 1; i <= a; ++i)
for (reg int j = i; j <= a; j += i)
fa[i] += a / j;
for (reg int i = 1; i <= b; ++i)
for (reg int j = i; j <= b; j += i)
fb[i] += b / j;
for (reg int i = 1; i <= c; ++i)
for (reg int j = i; j <= c; j += i)
fc[i] += c / j;
// u = v = w
LL ans = 0;
for (reg int i = 1; i <= mn; ++i)
if (mu[i])
ans += mu[i] * mu[i] * mu[i] * fa[i] * fb[i] * fc[i];
// u, v, w 有两个相等
etot = 0;
for (reg int i = 1; i <= mx; i += 4) {
deg[i] = 0; deg[i + 1] = 0; deg[i + 2] = 0; deg[i + 3] = 0;
}
for (reg int g = 1; g <= mx; ++g)
for (reg int i = 1; i * g <= mx; ++i)
if (mu[i * g])
for (reg int j = i + 1; 1ll * i * j * g <= mx; ++j)
if (mu[j * g] && gcd(i, j) == 1) {
int u = i * g, v = j * g, lcm = i * j * g;
int uuv = mu[u] * mu[u] * mu[v], uvv = mu[u] * mu[v] * mu[v];
ans += uuv * (fa[u] * fb[lcm] * fc[lcm] + fa[lcm] * fb[u] * fc[lcm] + fa[lcm] * fb[lcm] * fc[u]);
ans += uvv * (fa[v] * fb[lcm] * fc[lcm] + fa[lcm] * fb[v] * fc[lcm] + fa[lcm] * fb[lcm] * fc[v]);
++deg[u]; ++deg[v];
es[++etot] = Edge(u, v, lcm);
}
// u, v, w 两两不同
for (reg int i = 1; i <= mx; ++i) G[i].clear();
for (reg int i = 1; i <= etot; ++i) {
int u = es[i].u, v = es[i].v, w = es[i].w;
if (deg[u] > deg[v] || (deg[u] == deg[v] && u < v))
G[u].push_back(Node(v, w));
else G[v].push_back(Node(u, w));
}
for (reg int x = 1; x <= mx; ++x) if (mu[x]) {
unsigned sz = G[x].size();
for (reg unsigned i = 0; i < sz; ++i)
vis[G[x][i].to] = G[x][i].w;
for (reg unsigned i = 0; i < sz; ++i) {
int y = G[x][i].to, wxy = G[x][i].w;
if (!mu[y]) continue;
for (reg unsigned j = 0, sszz = G[y].size(); j < sszz; ++j) {
int z = G[y][j].to;
if (vis[z]) {
int wyz = G[y][j].w, wxz = vis[z], mmuu = mu[x] * mu[y] * mu[z];
if (!mmuu) continue;
ans += mmuu * (fa[wxy] * fb[wyz] * fc[wxz]
+ fa[wxy] * fb[wxz] * fc[wyz]
+ fa[wyz] * fb[wxy] * fc[wxz]
+ fa[wyz] * fb[wxz] * fc[wxy]
+ fa[wxz] * fb[wxy] * fc[wyz]
+ fa[wxz] * fb[wyz] * fc[wxy]);
}
}
}
for (reg unsigned i = 0; i < sz; ++i)
vis[G[x][i].to] = 0;
}
writeln(ans % mod);
}
return 0;
}
另附 的做法:
设函数 那么
枚举 ,枚举 ,这个 直接 暴力就行了,因为 所以枚举 和计算 的总复杂度是 (调和级数乘以 ),总复杂度