My Blogs
杜教筛是一种能在 $\mathcal O(n^{\frac 2 3})$ 的时间复杂度内求积性函数前缀和的筛法。虽然复杂度比较优秀,但是被筛的积性函数需要满足特殊性质。
Min_25 筛由 Min_25 发明,相对更通用,其时间复杂度为 $\mathcal O(\frac{n^{\frac 3 4}}{\log n})$。
首先构造一个完全积性函数 $f’(x)$ 满足其在质数处和 $f(x)$ 取值相同,其余位置任意。下设 $\mathbb{P}$ 为质数集合,$p_j$ 为第 $j$ 个质数。
设 $g(n,j)$ 表示 $\sum_{i\leq n}f’(i)[i\in \mathbb{P}\lor\min_{x\vert i,x\in\mathbb{P}}>p_j]$。一开始 $g(n,0)$ 就是 $\sum_{i\leq n}f’(i)$,考虑 $j$ 变大了之后会减少那些数的贡献:
若 $p_j^2>n$,则显然 $g(n,j)=g(n,j-1)$。
否则,有贡献的是最小质因子恰为 $p_j$ 的数。可以考虑整体除掉一个 $p_j$,因为满足完全积性,所以所有应当删去的树贡献就全部少乘了 $f’(p_j)$。
这样就变成了 $g(\lfloor\frac n {p_j}\rfloor,j-1)$。但是这样回多算 $[1,\lfloor\frac n {p_j}\rfloor]$ 内质数的贡献,需要再减去。
转移方程为:
\[g(n,j)=\begin{cases} g(n,j-1)\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad\qquad p_j^2>n\\ g(n,j-1)-f'(p_j)(g(\lfloor\frac n {p_j}\rfloor,j-1)-\sum_{j\leq \lfloor\frac n {p_j}\rfloor,j\in\mathbb{P}}f'(j))\;\;\;\;\;\;p_j^2\leq n \end{cases}\]这样最后得到的 $g(n,\lvert\mathbb P\rvert)$ 就是 $n$ 以内所有质数处的函数值。经过上述操作可以把错误的函数值舍去,最后剩下的是对的函数值。
预处理 $g$ 之后设 $S(n,j)$ 表示 $\sum_{i\leq n}f(i)(\min_{p\vert i,p\in \mathbb P}p>p_j)$。
首先对于 $i$ 是质数的部分,答案就是 $g(n,\lvert\mathbb P\rvert)$,这样会多算 $<j$ 的质数,再减去 $\sum_{i\in\mathbb P,i<j} f(i)$。
否则枚举其最小的因子 $p$ 和次幂 $e$,则有:
\[S(n,j)=g(n,\lvert\mathbb P\rvert)-\sum_{i\in\mathbb P,i<j} f(i)+\sum_{i>j}\sum_{e\geq 1}^{p_i^e\leq n}f(p_i^e)(S(\lfloor\frac n {p_i^e}\rfloor,i)+[e>1])\]加入 $[e>1]$ 是因为质数处的点值已经被计算过,但是 $p^k(k>1)$ 的点值还没有被计算。
边界条件是 $S(n,j)=0(p_j>n)$。但是这样预处理 $g$ 的复杂度达到了 $\mathcal O(n\log n)$。
注意到 $g(i,j)$ 有用当且仅当存在 $k$ 使得 $\lfloor\frac n k\rfloor=i$,所以有用的点只有 $\mathcal O(\sqrt n)$ 个。总复杂度是 $\mathcal O(\frac{n^{\frac 3 4}}{\log n})$ 不会证。
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
ll n,w[300010],g1[300010],g2[300010],sum1[300010],sum2[300010],pr[300010];
int N,m,cnt,id1[300010],id2[300010],iv6;
bool v[300010];
inline void sieve()
{
for(int i=2;i<=N;++i)
{
if(!v[i])pr[++cnt]=i;
for(int j=1;j<=cnt&&pr[j]*i<=N;++j)
{v[pr[j]*i]=1;if(i%pr[j]==0)break;}
}
for(int i=1;i<=cnt;++i)
{
sum1[i]=Cadd(sum1[i-1],pr[i]);
sum2[i]=(sum2[i-1]+1ll*pr[i]*pr[i])%MOD;
}
}
inline int f1(ll x){x%=MOD;return x*(x+1)/2%MOD;}
inline int f2(ll x){x%=MOD;return x*(x+1)%MOD*(2*x+1)%MOD*iv6%MOD;}
inline int getid(ll x){return x<=N?id1[x]:id2[n/x];}
ll S(ll x,int j)
{
if(pr[j]>x)return 0;
ll ans=Cdel(Cdel(g2[getid(x)],g1[getid(x)]),Cdel(sum2[j],sum1[j]));
for(int i=j+1;i<=cnt&&pr[i]*pr[j]<=x;++i)
{
for(ll e=1,sp=pr[i];sp<=x;sp*=pr[i],++e)
ans=(ans+sp%MOD*(sp%MOD-1)%MOD*(S(x/sp,i)+(e>1)))%MOD;
}
return ans;
}
inline void mian()
{
read(n),N=sqrt(n),sieve(),iv6=power(6,MOD-2);
for(ll l=1,r;l<=n;l=r+1)
{
r=n/(n/l),w[++m]=n/l;
g1[m]=f1(w[m])-1,g2[m]=f2(w[m])-1;
if(w[m]<=N)id1[w[m]]=m;else id2[n/w[m]]=m;
}
for(int i=1;i<=cnt;++i)
{
for(int j=1;j<=m&&pr[i]*pr[i]<=w[j];++j)
{
g1[j]=((g1[j]-1ll*pr[i]*(g1[getid(w[j]/pr[i])]-sum1[i-1]))%MOD+MOD)%MOD;
g2[j]=((g2[j]-1ll*pr[i]*pr[i]%MOD*(g2[getid(w[j]/pr[i])]-sum2[i-1]))%MOD+MOD)%MOD;
}
}
cout<<S(n,0)+1;
}