P3327 [SDOI2015]约数个数和 莫比乌斯反演
P3327 [SDOI2015]约数个数和 莫比乌斯反演
链接
思路
第一个式子我也不会,luogu有个证明,自己感悟吧。
\[d(ij)=\sum\limits_{x|i}\sum\limits_{y|j}[gcd(x,y)==1]
\]
\[\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}\sum\limits_{x|i}\sum\limits_{y|j}[gcd(x,y)==1]
\]
\[\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}{\left \lfloor \frac{n}{i} \right \rfloor \left \lfloor \frac{m}{j} \right \rfloor \left [ gcd(i,j)==1 \right ]}
\]
\[f(x)=\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}{\left \lfloor \frac{n}{i} \right \rfloor \left \lfloor \frac{m}{j} \right \rfloor \left [ gcd(i,j)==x \right ]}
\]
\[g(x)=\sum\limits_{x|d} f(d)
\]
\[g(x)=\sum\limits_{i=1}^{n}\sum\limits_{j=1}^{m}{\left \lfloor \frac{n}{i} \right \rfloor \left \lfloor \frac{m}{j} \right \rfloor \left [ x|gcd(i,j)\right ]}
\]
\[g(x)=\sum\limits_{i=1}^{\frac{n}{x}}\sum\limits_{j=1}^{\frac{m}{x}}{\left \lfloor \frac{n}{x*i} \right \rfloor \left \lfloor \frac{m}{x*j} \right \rfloor }
\]
\[g(x)=\sum\limits_{i=1}^{\frac{n}{x}}{\left \lfloor \frac{n}{x*i} \right \rfloor }\sum\limits_{j=1}^{\frac{m}{x}}{\left \lfloor \frac{m}{x*j} \right \rfloor }
\]
\[g(x)=\sum\limits_{i=1}^{N}{\left \lfloor \frac{N}{i} \right \rfloor }\sum\limits_{j=1}^{M}{\left \lfloor \frac{M}{j} \right \rfloor }(N=n/x,M=m/x)
\]
整除分块预处理,O(1)查询g(x)
\[f(x)=\sum\limits_{x|d}\mu(\frac{d}{n})g(d)
\]
所求$$f(1)=\sum\limits_{d=1}^{min(m,n)}\mu(d)g(d)$$
g是可以整除分块的
其他
改马蜂,加空格
代码
#include <bits/stdc++.h>
const int N = 5e5+7;
using namespace std;
int read() {
int x = 0, f = 1; char s = getchar();
for (;s > '9' || s < '0'; s = getchar()) if (s == '-') f = -1;
for (;s >= '0' && s <= '9'; s = getchar()) x = x * 10 + s - '0';
return x * f;
}
int n, m, T;
int pri[N], vis[N], tot, mu[N], g[N];
void Euler(int limit) {
mu[1] = 1;
for (int i = 2; i <= limit; ++i) {
if (!vis[i]) {
pri[++tot] = i;
mu[i] = -1;
}
for (int j = 1; j <= tot && i * pri[j] <= limit; ++j) {
vis[i * pri[j]] = 1;
if (i % pri[j] == 0) {
mu[i * pri[j]] = 0;
break;
}
mu[i * pri[j]] = -mu[i];
}
}
for (int i = 1; i <= limit; ++i) {
for (int l = 1, r; l <= i; l = r + 1) {
r = i / (i / l);
g[i] += (r - l + 1) * (i / r);
}
mu[i] += mu[i - 1];
}
}
void solve() {
n = read(), m = read();
if (n > m) swap(n, m);
long long ans = 0;
for (int l = 1, r; l <= n; l = r + 1) {
r = min(n / (n / l), m / (m / l));
ans += 1LL * (mu[r] - mu[l-1]) * (1LL * g[n/l] * g[m/l]);
}
printf("%lld\n", ans);
}
int main() {
Euler(50000);
int T = read();
while (T--) solve();
return 0;
}