杜教筛

杜教筛

杜教筛用来解决类似求$\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)$时需要记忆化。记忆化有两种实现方式:

  1. 使用unordered_map
  2. 由于$\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$的前缀和

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
#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;
}
Your browser is out-of-date!

Update your browser to view this website correctly. Update my browser now

×