杜教筛
最近在透彻一些数学知识,多项式啊,线性代数啊,杜教筛也逃不掉。。。
洛咕博客的题解都会搬过来的,然后我会整理一遍所有文章。。。
然后下面是正题
我们现在已知积性函数\(f(n)\),求\(f(n)\)的前n项和\(S(n)\)
由于\(f(n)\)是积性函数,所以可以通过线性筛在\(O(n)\)的复杂度解决
但是由于有杜教筛这种东西使得复杂度被进一步优化
我们先扯扯Dirichlet卷积,假设\(f(n)\)和\(g(n)\)是两个数论函数,则\(f(n)\)和\(g(n)\)的Dirichlet卷积定义为\(\displaystyle(f*g)(n)=\sum_{d|n}f(d)g\left(\frac{n}{d}\right)\)或者写一个直观点的\(\displaystyle (f*g)(n)=\sum_{ab=n}f(a)g(b)\)
然后呢有一个性质就是如果f和g都是积性函数,他们的Dirichlet卷积也是积性函数。。。
然后呢我们假设\(h(n)=(f*g)(n)\)然后呢推式子
先是\(h(n)\)的前缀和式子
\(\displaystyle \sum_{i=1}^nh(i)=\sum_{i=1}^ng(i)S\left(\left\lfloor\frac n i\right\rfloor\right)\)
这是什么鬼woc
\(\displaystyle\sum_{i=1}^nh(i)\)
写成Dirichlet卷积的形式
\(\displaystyle\sum_{i=1}^n\sum_{d|i}g(d)f\left(\frac i d\right)\)
把d提取出来
\(\displaystyle \sum_{d=1}^ng(d)\sum_{d|i}f\left(\frac i d \right)\)
i是d的倍数,考虑让i除以d
\(\displaystyle \sum_{d=1}^ng(d)\sum_{i=1}^{\left\lfloor\frac n d \right \rfloor}f( i)\)
然后呢吧后面表示为前缀和的性质(把d字母再替换一下)
\(\displaystyle \sum_{i=1}^ng(i)S\left(\left\lfloor\frac n i\right\rfloor\right)\)
然后呢既然推出了这\(h(n)\)前缀和的另一种表示形式
那么我们先写一个无聊的式子\(\displaystyle g(1)S(n)=\sum_{i=1}^ng(i)S\left(\left\lfloor\frac n i\right\rfloor\right)-\sum_{i=2}^ng(i)S\left(\left\lfloor\frac n i\right\rfloor\right)\)
你会问我们推出来这么无聊的东西有什么用,把后面移项不就成恒等式了么
别急,右面的第一项不是上面推得\(h(i)\)的前缀和么。。。
所以\(\displaystyle g(1)S(n)=\sum_{i=1}^nh(i)-\sum_{i=2}^ng(i)S\left(\left\lfloor\frac n i\right\rfloor\right)\)
观察这个式子,就是\(S\)的递推式喽
根据数论分块思想 递归求S
然后呢...注意每求出一个前缀和用map存一下
一般我们做杜教筛时候,要自己钦定一个\(g\)和\(h=f*g\),使得\(g(x)\)能在O(1)获得,\(h(x)\)的前缀和能在\(O(1)\)获得,其实一般\(g(x)=1\)
然后呢求欧拉函数\(\varphi\)的前缀和,我们令\(g(x)=1\),则\(h(x)=x\)也就是id
所以\(\displaystyle S(n)=\sum_{i=1}^ni-\sum_{i=2}^nS\left(\left\lfloor\frac n i\right\rfloor\right)\)
然后呢根据等差数列有关知识\(\displaystyle \sum_{i=1}^ni=\frac{n(n+1)}2\),可以\(O(1)\)求出
后面就直接递归就行了
然后呢是莫比乌斯函数\(\mu\),其实更简单,我们还是令\(g(x)=1\),则\(h(x)=[x=1]\),也就是Dirichlet卷积得单位元函数
所以\(\displaystyle S(n)=1-\sum_{i=2}^nS\left(\left\lfloor\frac n i\right\rfloor\right)\)
然后呢说实现
这个算法时间复杂度是\(O(n^\frac{3}{4})\)
然后呢我不会证
然后据说先线性筛出前面的一些值可以把复杂度优化为\(O(n^{\frac{2}{3}})\)
然后呢这是一份洛谷板子题的代码
#include <bits/stdc++.h>
#define maxn 3000010
using namespace std;
map<int, long long> ans_phi;
map<int, int> ans_mu;
bool vis[maxn];
int prime[maxn], tot;
long long phi[maxn];
int mu[maxn];
void prework()
{
phi[1] = 1;
mu[1] = 1;
for (int i = 2; i <= 3000000; i++)
{
if (vis[i] == 0)
{
prime[++tot] = i;
mu[i] = -1;
phi[i] = i - 1;
}
for (int j = 1; j <= tot && prime[j] * i <= 3000000; j++)
{
vis[i * prime[j]] = 1;
if (i % prime[j] == 0)
{
mu[i * prime[j]] = 0;
phi[i * prime[j]] = phi[i] * prime[j];
break;
}
else
{
mu[i * prime[j]] = -mu[i];
phi[i * prime[j]] = phi[i] * (prime[j] - 1);
}
}
mu[i] += mu[i - 1];
phi[i] += phi[i - 1];
}
}
long long S_phi(int n)
{
if (n <= 3000000)
return phi[n];
if (ans_phi.count(n))
return ans_phi[n];
long long ans = 1LL * (n + 1) * n / 2;
for (int l = 2, r; l <= n; l = r + 1)
{
r= n / (n / l);
ans -= (r - l + 1) * S_phi(n / l);
}
return ans_phi[n] = ans;
}
int S_mu(int n)
{
if (n <= 3000000)
return mu[n];
if (ans_mu.count(n))
return ans_mu[n];
int ans = 1;
for (int l = 2, r; l <= n; l = r + 1)
{
r = n / (n / l);
ans -= (r - l + 1) * S_mu(n / l);
}
return ans_mu[n] = ans;
}
void read(int &x)
{
static char ch;
x = 0;
ch = getchar();
while (!isdigit(ch))
ch = getchar();
while (isdigit(ch))
{
x = x * 10 + ch - 48;
ch = getchar();
}
}
int main()
{
prework();
int T;
read(T);
for (int n, i = 1; i <= T; i++)
{
read(n);
printf("%lld %d\n", S_phi(n), S_mu(n));
}
return 0;
}
参考资料:洛谷题解