【luogu 1457】在表格里造序列（莫反）（杜教筛）

在表格里造序列

题目链接：luogu 1457

题目大意

有一个 n*n 的表格，每个位置有一些数列，其中 i 行 j 列里面的数列是所有满足长度为 m，值域是 1~n，且所有值的 gcd 与 i,j 的 gcd 相同的数组。
问你表格里总共有多少个数组。

思路

\(\color{white}{害，世界上最好的还得是莫反。}\)
\(\color{white}{当其它所有的题都在用奇奇怪怪的方式恶心你的时候}\)
\(\color{white}{只有莫反，纯纯的推式子，纯纯的，多好啊。}\)
\(\color{white}{（}\)
\(\color{white}{还有捏嘛csp怎么不考数学（算了考了我也不会，还是别考了吧}\)

---

考虑枚举 \(\gcd\) 的值，再分别枚举由多少个数组以及多少个表格上的位置满足条件。
\(\sum\limits_{d=1}^n(\sum\limits_{x=1}^n\sum\limits_{y=1}^n[\gcd(x,y)=d])(\sum\limits_{a_1=1}^n...\sum\limits_{a_m=1}^n[\gcd(a_1,...,a_m=d)])\)
\(\sum\limits_{d=1}^n(\sum\limits_{x=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum\limits_{y=1}^{\left\lfloor\frac{n}{d}\right\rfloor}[\gcd(x,y)=1])(\sum\limits_{a_1=1}^{\left\lfloor\frac{n}{d}\right\rfloor}...\sum\limits_{a_m=1}^{\left\lfloor\frac{n}{d}\right\rfloor}[\gcd(a_1,...,a_m=1)])\)
\(\sum\limits_{d=1}^n(\sum\limits_{x=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum\limits_{y=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum\limits_{p|\gcd(x,y)}\mu(p))(\sum\limits_{a_1=1}^{\left\lfloor\frac{n}{d}\right\rfloor}...\sum\limits_{a_m=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\sum\limits_{p|\gcd(a_1,...,a_m)}\mu(p))\)
\(\sum\limits_{d=1}^n(\sum\limits_{p=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\mu(p)\sum\limits_{x=1}^{\left\lfloor\frac{n}{dp}\right\rfloor}\sum\limits_{y=1}^{\left\lfloor\frac{n}{dp}\right\rfloor})(\sum\limits_{p=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\mu(p)\sum\limits_{a_1=1}^{\left\lfloor\frac{n}{dp}\right\rfloor}...\sum\limits_{a_m=1}^{\left\lfloor\frac{n}{dp}\right\rfloor})\)
\(\sum\limits_{d=1}^n(\sum\limits_{p=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\mu(p)(\left\lfloor\frac{n}{dp}\right\rfloor)^2)(\sum\limits_{p=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\mu(p)(\left\lfloor\frac{n}{dp}\right\rfloor)^m)\)

\(\sum\limits_{d=1}^n\sum\limits_{p=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\mu(p)(\left\lfloor\frac{n}{dp}\right\rfloor)^m=n^m\)
\(\sum\limits_{p=1}^n\mu(p)(\left\lfloor\frac{n}{p}\right\rfloor)^m=n^m-\sum\limits_{d=2}^n\sum\limits_{p=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\mu(p)(\left\lfloor\frac{n}{dp}\right\rfloor)^m\)
设 \(\sum\limits_{p=1}^n\mu(p)(\left\lfloor\frac{n}{p}\right\rfloor)^m=f_n\)

\(f_n=n^m-\sum\limits_{d=2}^nf_{\left\lfloor\frac{n}{d}\right\rfloor}\)

\(\sum\limits_{d=1}^n(\sum\limits_{p=1}^{\left\lfloor\frac{n}{d}\right\rfloor}\mu(p)(\left\lfloor\frac{n}{dp}\right\rfloor)^2)(f_{\left\lfloor\frac{n}{d}\right\rfloor})\)

其中 \(f\) 可以记忆化一下在常规整除分块做到 \(O(n^{\frac{3}{4}})\)，左边的一坨式子，都可以用杜教筛快速 rush \(\mu\) 前缀和也是 \(O(n^{\frac{3}{4}})\) 从而通过题目捏。

代码

#include<cmath>
#include<cstdio>
#define ll long long
#define mo 998244353

using namespace std;

const int Lim = 1e6;
const int B = 2e5 + 100;
ll n, m, lim, mus[Lim + 100];
ll rem_f[B], rem_mu[B], prime[Lim + 100];
bool np[Lim + 100];

int id(ll x) {
	if (x <= lim) return x; return lim + n / x;
}

ll ksm(ll x, ll y) {
	ll re = 1;
	while (y) {
		if (y & 1) re = re * x % mo;
		x = x * x % mo; y >>= 1;
	}
	return re;
}

ll get_f(ll now) {
	if (rem_f[id(now)] != -1) return rem_f[id(now)];
	ll re = ksm(now, m);
	for (ll L = 2, R; L <= now; L = R + 1) {
		R = now / (now / L);
		(re += mo - (R - L + 1) % mo * get_f(now / L) % mo) %= mo;
	}
	return rem_f[id(now)] = re;
}

ll ask(ll now) {
	if (now <= Lim) return (mus[now] % mo + mo) % mo;
	if (rem_mu[id(now)] != -1) return rem_mu[id(now)];
	ll re = 1;
	for (int L = 2, R; L <= now; L = R + 1) {
		R = now / (now / L);
		(re += mo - (R - L + 1) * ask(now / L) % mo) %= mo;
	}
	return rem_mu[id(now)] = re;
}

ll slove(ll n) {
	ll ans = 0;
	for (ll L = 1, R; L <= n; L = R + 1) {
		R = n / (n / L);
		ll mus = (ask(R) - ask(L - 1) + mo) % mo;
		(ans += mus * (n / L) % mo * (n / L) % mo) %= mo;
	}
	return ans;
}

int main() {
	for (int i = 0; i < B; i++) rem_f[i] = rem_mu[i] = -1;
	mus[1] = 1;
	for (int i = 2; i <= Lim; i++) {
		if (!np[i]) mus[i] = -1, prime[++prime[0]] = i;
		for (int j = 1; j <= prime[0] && i * prime[j] <= Lim; j++) {
			np[i * prime[j]] = 1;
			if (i % prime[j] == 0) break;
				else mus[i * prime[j]] = -mus[i];
		}
		mus[i] += mus[i - 1];
	}
	
	scanf("%lld %lld", &n, &m); lim = sqrt(n);
	
	ll ans = 0;
	for (ll L = 1, R; L <= n; L = R + 1) {
		R = n / (n / L);
		(ans += (R - L + 1) % mo * slove(n / L) % mo * get_f(n / L) % mo) %= mo;
	}
	printf("%lld", ans);
	
	return 0;
}