2020-07-25 5 分钟读完 (大约 744 个字)

杜教筛

杜教筛用来解决类似求$\sum_{i=1}^nf(i)$的问题

设$S(n)=\sum_{i=1}^nf(i)$

假如我们现在找到了一个积性函数$g$，那么考虑下式（其中$*$表示狄利克雷卷积）

$$
\sum_{i=1}^n(f \ast g)(i) \
=\sum_{i=1}^n\sum_{d|i}f(\frac{i}{d})g(d) \
=\sum_{d=1}^ng(d)\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}f(i) \
=\sum_{d=1}^ng(d)S(\lfloor\frac{n}{d}\rfloor)
$$

再考虑下式

$$
g(1)S(n)=\sum_{i=1}^ng(i)S(\lfloor\frac{n}{i}\rfloor)-\sum_{i=2}^ng(i)S(\lfloor\frac{n}{i}\rfloor) \
=\sum_{i=1}^n(f \ast g)(i)-\sum_{i=2}^ng(i)S(\lfloor\frac{n}{i}\rfloor)
$$

那么，假如我们可以快速地求出$(f*g)(i)$的前缀和以及$g(i)$的前缀和，那么我们就可以用数论分块递归地求解$S(n)$

其中，计算$S(n)$时需要记忆化。记忆化有两种实现方式：

使用unordered_map
由于$\lfloor\frac{\frac{n}{i}}{j}\rfloor \in {\lfloor\frac{n}{i}\rfloor}$，所以不同的$S(k)$的$k$的个数最多只有$O(\sqrt{n})$个。因此当$k\leq \sqrt{n}$时，记录在$sum(k)$中，否则记录在$sum(\sqrt{n}+\lfloor\frac{n}{k}\rfloor)$中

时间复杂度

假如直接数论分块递归求解，那么时间复杂度
$$
T(n)=\sum_{i=1}^{\sqrt{n}}\sqrt{i}+\sum_{i=1}^{\sqrt{n}}\sqrt{\frac{n}{i}}=O(n^{\frac{3}{4}})
$$
用积分等方法可以轻松计算

假如先预处理出前$m$个$S(n)$，那么时间复杂度变为
$$
T(n)=m+\sum_{i=1}^{\lfloor\frac{n}{m}\rfloor}\sqrt{\frac{n}{i}}=O(m+\frac{n}{\sqrt{m}})
$$

取$m=n^{\frac{2}{3}}$，即得
$$
T(n)=O(n^{\frac{2}{3}})
$$

$g$的选择

可以发现能否使用杜教筛的关键就在于$g$的选择，下列几种常见的$f$与$g$的对应

$f=\mu$

$$
f=\mu,g=I,f*g=\epsilon
$$

$f=\phi$

$$
f=\phi,g=I,f*g=id
$$

$f=\phi \cdot id$

$$
f=\phi \cdot id,g=id,f*g=id_2
$$

代码示例

求$\phi$和$\mu$的前缀和

#include <cstdio>
#include <unordered_map>

long long read()
{
    char last = '+', ch = getchar();
    while (ch < '0' || ch > '9') last = ch, ch = getchar();
    long long tmp = 0;
    while (ch >= '0' && ch <= '9') tmp = tmp * 10 + ch - 48, ch = getchar();
    if (last == '-') tmp = -tmp;
    return tmp;
}

const int lim = 5000000, _l = lim + 10;
int T;
int flag[_l], phi[_l], mu[_l];
long long sphi[_l], smu[_l];
int primenum;
int prime[_l];
std::unordered_map<int, long long> mphi, mmu;

inline long long Sphi(int n)
{
	if (n <= lim)
	{
		return sphi[n];
	}
	if (mphi[n])
	{
		return mphi[n];
	}
	long long ans = (long long)n * (n + 1) / 2;
	for (int l = 2, r; l <= n; l = r + 1)
	{
		r = n / (n / l);
		ans -= (r - l + 1) * Sphi(n / l);
	}
	return mphi[n] = ans;
}

inline long long Smu(int n)
{
	if (n <= lim)
	{
		return smu[n];
	}
	if (mmu[n])
	{
		return mmu[n];
	}
	long long ans = 1;
	for (int l = 2, r; l <= n; l = r + 1)
	{
		r = n / (n / l);
		ans -= (r - l + 1) * Smu(n / l);
	}
	return mmu[n] = ans;
}

int main()
{
	//init
	mu[1] = 1;
	phi[1] = 1;
	for (int i = 2; i <= lim; i++)
	{
		if (! flag[i])
		{
			primenum++;
			prime[primenum] = i;
			mu[i] = -1;
			phi[i] = i - 1;
		}
		for (int j = 1; j <= primenum && i * prime[j] <= lim; j++)
		{
			flag[i * prime[j]] = 1;
			if (i % prime[j] == 0)
			{
				mu[i * prime[j]] = 0;
				phi[i * prime[j]] = phi[i] * prime[j];
			}
			else
			{
				mu[i * prime[j]] = -mu[i];
				phi[i * prime[j]] = phi[i] * (prime[j] - 1);
			}
		}
	}
	for (int i = 1; i <= lim; i++)
	{
		smu[i] = smu[i - 1] + mu[i];
		sphi[i] = sphi[i - 1] + phi[i];
	}

	int T = read();
	while (T--)
	{
		int n = read();
		printf("%lld %lld\n", Sphi(n), Smu(n));
	}
	return 0;
}

# 数学, 数论, 算法

杜教筛

杜教筛