[牛客挑战赛46E] 反演（莫比乌斯反演）

题意

分析

首先，根据 ${\sigma }_0(ij)=\sum _{x|i}\sum _{y|j}[(x,y)==1]$，有
$$\begin{gathered} F(m,n) \\ =\sum _{i=1}^n\sum _{j|m}\sum _{x|i}\sum _{y|j}[(x,y)==1] \\ =\sum _{x=1}^n\sum _{y|m}\lfloor \frac{n}{x}\rfloor {\sigma }_0(\frac{m}{y})[(x,y)==1] \\ =\sum _{x=1}^n\sum _{y|m}\lfloor \frac{n}{x}\rfloor {\sigma }_0(\frac{m}{y})\sum _{d|x,d|y}\mu (d) \\ =\sum _{d|m}\mu (d)h(\lfloor \frac{n}{d}\rfloor )f(\frac{m}{d}) \end{gathered}$$
其中，$h(n)=\sum _{i=1}^n\lfloor \frac{n}{i}\rfloor$，$f(n)=\sum _{d|n}{\sigma }_0(d)$。
那么，我们要求的是 $F(m!,n)$，即 $\sum _{d|m!}\mu (d)h(\lfloor \frac{n}{d}\rfloor )f(\frac{m!}{d})$。
可以考虑枚举 $d$，因为 $100$ 以内有 $25$ 个质数，因此最多有 $2^{25}$ 个 $d$。
因为 $f$ 是个积性函数，可以快速求出。
每得到一个 $d$，除法分块得到 $h$，用哈希表储存值，这部分复杂度是 $O(n^{\frac{3}{4}})$。
总的复杂度是 $O(2^{25}+n^{\frac{3}{4}})$。

代码

#include <bits/stdc++.h>
#include <ext/pb_ds/assoc_container.hpp>
#include <ext/pb_ds/hash_policy.hpp>
using namespace std;
using namespace __gnu_pbds;
typedef long long LL;

void debug_out(){
    cerr << endl;
}
template<typename Head, typename... Tail>
void debug_out(Head H, Tail... T){
    cerr << " " << to_string(H);
    debug_out(T...);
}
#ifdef local
#define debug(...) cerr<<"["<<#__VA_ARGS__<<"]:",debug_out(__VA_ARGS__)
#else
#define debug(...) 55
#endif

const int mod = 998244353, N = 1e7 + 5;
LL n;
int p[26], cnt, power[105], inv[N], ans, tot = 1;
gp_hash_table<int, int> f;
int h(LL n){
	if(f[n]) return f[n];
	LL s = 0;
	for(LL l = 1, r; l <= n; l = r + 1){
		r = n / (n / l);
		s += (r - l + 1) * (n / l) % mod;
		if(s >= mod) s -= mod;
	}
	return f[n] = s;
}
void dfs(int pos, LL x, int mu, int sum){
	if(x > n) return;
	if(pos == cnt + 1){
		ans = (ans + (LL)mu * h(n / x) * sum % mod) % mod;
		return;
	}
	dfs(pos + 1, x, mu, sum);
	dfs(pos + 1, x * p[pos], -mu, (LL)sum * inv[(power[p[pos]] + 2) * (power[p[pos]] + 1) / 2] % mod * ((power[p[pos]] + 1) * power[p[pos]] / 2) % mod);
}
int main(){
	ios::sync_with_stdio(false);
	cin.tie(0), cout.tie(0);
	int m;
	inv[1] = 1;
	for(int i = 2; i <= 10000; i++) inv[i] = (LL)inv[mod % i] * (mod - mod / i) % mod;
	cin >> n >> m;
	for(int i = 2; i <= m; i++){
		int ok = 1;
		for(int j = 2; j <= i / j; j++) if(i % j == 0) ok = 0;
		if(ok){
			p[++cnt] = i;
			int t = m;
			while(t){
				power[i] += t / i;
				t /= i;
			}
			tot = (LL)tot * ((power[i] + 2) * (power[i] + 1) / 2) % mod;
		}
	}
	dfs(1, 1, 1, tot);
	cout << (ans + mod) % mod << '\n';
	return 0;
}