min25筛

min_25 筛

非常棒的博客

已知一积性函数 \(f(x)\),\(f(p)\) 为关于 \(p\) 的多项式,且 \(f(p^k)\) 可以快速算出来,求 \(\sum_{i=1}^nf(i)\)。

例:求 \(\sum_{i=1}^n \sigma_0(i^k)\),其中 \(n = 10^{10}\)

step1 : 求质数的积性函数前缀和

求 \(1 ... n\) 以内的质数个数。

我们设 \(f(n,m)\) 表示 \(1...n\) 中的所有质数或者不含最小质因子为 \(p_1...p_m\) (\(p_i\) 表示第 \(i\) 个质数)的数的个数。

那么有:

\[f(n,m) = f(n,m-1) + (f(\lfloor \frac{n}{p_m} \rfloor,m-1) - (m-1)) \]

意义为:考虑从 \(f(n,m-1)\) 处转移过来,那么我们只需要筛掉最小质因子恰好是 \(p_m\) 的数。那么我们把所有要筛的数提出来个 \(p_m\),剩下的数的最小质因子为 \(p_1...p_{m-1}\)。注意删掉“质数"部分。

更普遍的情况。

类似地,我们设 \(f(n,m)\) 表示 \(1...n\) 中的所有质数或者不含最小质因子为 \(p_1...p_m\) (\(p_i\) 表示第 \(i\) 个质数)的数的函数和。这时候我们需要认为合数也有一个假的和质数一样的函数和。那么有:

\[f(n,m) = f(n,m-1) + F(p_m)(f(\lfloor \frac{n}{p_m} \rfloor,m-1) - psum(m-1)) \]

(函数重名了,用 \(F(n)\) 代替了积性函数;\(psum\) 为预处理出来的质数函数前缀和)

由于是多项式,我们可以拆成 \(F(p) = p^k\) 的形式,最终拼起来。这样,假函数成为完全积性函数,直接提出来 \(F(p_m)\) 就是对的了。

复杂度证明

见那篇博客。

发现 \(p_m^2 \le n\)

可以积分证明,复杂度为 \(O(n^{0.75} / \log n)\)

step2 : 求所有数的积性函数和

我们设 \(g(n,m)\) 表示 \(1...n\) 中最小质因子大于 \(p_m\) 的数的(真)积性函数和。

那么有:

\[g(n,m) = f(n,m) - psum(m-1) + \sum_{p_k^e} F(p_k^e)(g(\lfloor \frac{n}{p^e} \rfloor,k) + [e \not= 1]) \]

意义为:分质数,合数考虑贡献。质数的贡献我们已经搞好了,为 \(f(n,m) - psum(m-1)\),合数的话我们需要枚举最小质因子及其次数,即枚举 \(k > m\) 的 \(p_k\),再枚举次数 \(e\)。由于是积性函数,我们可以直接将 \(F(p_k^e)\) 提出来。值得注意的是,形如 \(p^k,k > 1\) 的数我们既没有在质数处枚举到,也没有在合数处提出来 \(p_k^e\) 后枚举到,所以需要加上个 \([e \not= 1]\)。

复杂度:不会证。状态数和上面一样是 \(O(n^{0.75} / \log n)\) 的,但是转移复杂度比较玄学,据说总复杂度是 \(O(n^{0.75} / \log n)\) 的,也有的说是 \(O(n^{1 - \epsilon})\) 的。

代码:

//divcntk
int block, itot, pos[N];
ull id[N], n, K, f[N];
inline int get_id(ull x) {
	return pos[x <= block ? x : block + n / x];
}
inline void get_f() {
	for (ull l = 1, r; l <= n; l = r + 1) {
		ull zhi = n / l;
		r = n / zhi;
		id[++itot] = zhi;
		pos[zhi <= block ? zhi : block + n / zhi] = itot;
		f[itot] = zhi - 1;
	}
	for (int j = 1; 1ull * pri[j] * pri[j] <= n; ++j) {
		for (int i = 1; id[i] >= 1ull * pri[j] * pri[j]; ++i) {
			f[i] = f[i] - f[get_id(id[i] / pri[j])] + j - 1;
		}
	}
}

ull get_g(ull n, int j) {
	if (pri[j] >= n)	return 0;
	ull res = (f[get_id(n)] - j) * (K + 1);//bug
	for (int k = j + 1; 1ull * pri[k] * pri[k] <= n; ++k) {
		ull tmp = pri[k];
		for (int e = 1; tmp <= n; ++e, tmp = tmp * pri[k]) {
			res += (K * e + 1) * (get_g(n / tmp, k) + (e > 1));
		}
	}
	return res;
}

inline void work() {
	read(n), read(K); block = sqrt(n);//bug
	get_f();
	ull res = get_g(n, 0) + 1;
	printf("%llu\n", res);
}

min_25 筛

上一篇:个人笔记·进制哈希


下一篇:豪华版点名器