【题目】
求积性函数
f
(
p
c
)
=
p
⨁
c
f(p^c)=p\bigoplus c
f(pc)=p⨁c的前缀和,其中
p
p
p是素数,
f
(
1
)
=
1
f(1)=1
f(1)=1
【解题思路】
学了一波
Min_25
\text{Min\_25}
Min_25筛。
这个东西可以在
O
(
n
3
4
log
n
)
O(\frac {n^{\frac 3 4} }{\log n})
O(lognn43)的复杂度内求出积性函数
f
(
x
)
f(x)
f(x)的前缀和,同时要求对于素数
p
p
p,
f
(
p
)
f(p)
f(p)是一个关于
p
p
p的简单多项式,
f
(
p
c
)
f(p^c)
f(pc)可以快速计算。
主要的思想是将计算拆成素数部分和非素数部分计算。
详细过程及推导可以看这里:blog
较为简单地说一下过程:
第一步要对每个
x
=
⌊
n
i
⌋
x=\lfloor \frac n i \rfloor
x=⌊in⌋,求出
∑
i
=
1
x
[
i
is prime
]
f
(
i
)
\sum_{i=1}^x [i\text{ is prime}]f(i)
∑i=1x[i is prime]f(i)。
先筛出
n
\sqrt n
n内所有素数,我们设
P
i
P_i
Pi表示从小到大第
i
i
i个素数。
设
g
(
n
,
j
)
=
∑
i
=
1
n
[
i
∈
P
o
r
min
(
p
)
>
P
j
]
f
(
i
)
g(n,j)=\sum_{i=1}^{n}[i \in P \ or\ \min(p)>P_j]f(i)
g(n,j)=∑i=1n[i∈P or min(p)>Pj]f(i)
即
i
i
i是质数,或者
i
i
i的最小质因子大于
P
j
P_j
Pj时,统计入
g
(
n
,
j
)
g(n,j)
g(n,j)
那么我们可以发现我们上面要求的哪个东西就是
g
(
x
,
∣
P
∣
)
g(x,|P|)
g(x,∣P∣),而
g
g
g有两种转移。
g
(
n
,
j
)
=
{
g
(
n
,
j
−
1
)
P
j
2
>
n
g
(
n
,
j
−
1
)
−
f
(
P
j
)
[
g
(
n
P
j
,
j
−
1
)
−
∑
i
=
1
j
−
1
f
(
P
i
)
]
P
j
2
≤
n
g(n,j)=\begin{cases} g(n,j-1)&P_j^2\gt n\\ g(n,j-1)-f(P_j)[g(\frac{n}{P_j},j-1)-\sum_{i=1}^{j-1}f(P_i)]&P_j^2\le n\end{cases}
g(n,j)={g(n,j−1)g(n,j−1)−f(Pj)[g(Pjn,j−1)−∑i=1j−1f(Pi)]Pj2>nPj2≤n
初始时
g
(
n
,
0
)
g(n,0)
g(n,0)表示将所有数当作素数代入
f
(
p
)
f(p)
f(p)的和。
接下来设
S
(
n
,
j
)
=
∑
i
=
1
n
[
min
(
p
)
≥
P
j
]
f
(
i
)
S(n,j)=\sum_{i=1}^n[\min(p)\ge P_j]f(i)
S(n,j)=∑i=1n[min(p)≥Pj]f(i),即所有最小质因子大于等于
P
j
P_j
Pj的
f
f
f之和。
最后的答案就是
S
(
n
,
1
)
+
f
(
1
)
S(n,1)+f(1)
S(n,1)+f(1)。
素数部分的答案我们已经有了,即
g
(
n
,
j
)
−
∑
i
=
1
j
−
1
f
(
P
i
)
g(n,j)-\sum_{i=1}^{j-1}f(P_i)
g(n,j)−∑i=1j−1f(Pi)。(保证最小质因子大于等于
P
j
P_j
Pj所以要把小于它的质数减掉)
考虑合数,枚举它最小质因子及其出现次数,由于
f
f
f是积性的,所以可以直接相乘。
S
(
n
,
j
)
=
g
(
n
,
j
)
−
∑
i
=
1
j
−
1
f
(
P
i
)
+
∑
k
=
j
P
k
2
≤
n
∑
e
=
1
P
k
e
+
1
≤
n
S
(
n
P
k
e
,
k
+
1
)
×
f
(
P
k
e
)
+
f
(
P
k
e
+
1
)
S(n,j)=g(n,j)-\sum_{i=1}^{j-1}f(P_i)+\sum_{k=j}^{P_k^2\le n}\sum_{e=1}^{P_k^{e+1}\le n}S(\frac{n}{P_k^e},k+1)\times f(P_k^e)+f(P_k^{e+1})
S(n,j)=g(n,j)−i=1∑j−1f(Pi)+k=j∑Pk2≤ne=1∑Pke+1≤nS(Pken,k+1)×f(Pke)+f(Pke+1)
对于这一道题,我们有积性函数
f
(
p
c
)
=
p
⨁
c
f(p^c)=p\bigoplus c
f(pc)=p⨁c,求前
n
n
n项和。
可以观察到除了
2
2
2以外的素数都满足
f
(
p
)
=
p
−
1
f(p)=p-1
f(p)=p−1,而
f
(
2
)
=
2
+
1
=
3
f(2)=2+1=3
f(2)=2+1=3,因此我们分别计算出
g
(
x
,
∣
P
∣
)
=
∑
i
=
1
x
[
i
is prime
]
i
g(x,|P|)=\sum_{i=1}^x[i\text{ is prime}]i
g(x,∣P∣)=∑i=1x[i is prime]i和
h
(
x
,
∣
P
∣
)
=
∑
i
=
1
x
[
i
is prime
]
1
h(x,|P|)=\sum_{i=1}^x[i\text{ is prime}]1
h(x,∣P∣)=∑i=1x[i is prime]1,那么可以得到
S
(
n
,
j
)
S(n,j)
S(n,j)。
S
(
n
,
j
)
=
g
(
n
,
∣
P
∣
)
−
∑
i
=
1
j
−
1
(
P
i
−
1
)
−
h
(
n
,
∣
P
∣
)
S(n,j)=g(n,|P|)-\sum_{i=1}^{j-1}(P_i-1)-h(n,|P|)
S(n,j)=g(n,∣P∣)−i=1∑j−1(Pi−1)−h(n,∣P∣)
这里首先是所有素数的
f
(
x
)
f(x)
f(x)的值,即
g
g
g函数,最小素数为
P
j
P_j
Pj,所以减去前面的贡献,然后
f
(
x
)
=
x
−
1
f(x)=x-1
f(x)=x−1,减去多算的那个
1
1
1。
注意当
j
=
1
j=1
j=1特判上
+
2
+2
+2。
xjb构造了一种权值对应关系避免了注释中的开两个数组,但是相对麻烦
【参考代码】
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N=1e6+10,mod=1e9+7;
int lim,m,pnum,pri[N],bo[N];
ll n,w[N],g[N],h[N],sq[N];
void initpri()
{
bo[1]=1;
for(int i=2;i<=lim;++i)
{
if(!bo[i]) pri[++pnum]=i,sq[pnum]=(sq[pnum-1]+i)%mod;
for(int j=1;j<=pnum && i*pri[j]<=lim;++j)
{
bo[i*pri[j]]=1;
if(!(i%pri[j])) break;
}
}
}
ll id(ll x){return x<=lim?m-x+1:(n/x);}
ll S(ll x,int y)
{
if(x<=1 || pri[y]>x) return 0;
//int k=(x<=lim)?id1[x]:id2[n/x];
int k=id(x);ll res=((g[k]-h[k]-sq[y-1]+y-1+(y==1?1:0)*2)%mod+mod)%mod;
for(int i=y;i<=pnum && (ll)pri[i]*pri[i]<=x;++i)
{
ll p1=pri[i],p2=(ll)pri[i]*pri[i];
for(int c=1;p2<=x;++c,p1=p2,p2*=pri[i])
res=(res+S(x/p1,i+1)*(pri[i]^c)%mod+(pri[i]^(c+1)))%mod;
}
return res;
}
int main()
{
#ifndef ONLINE_JUDGE
freopen("LOJ6053.in","r",stdin);
freopen("LOJ6053.out","w",stdout);
#endif
scanf("%lld",&n);lim=sqrt(n);initpri();
for(ll i=1,j;i<=n;i=j+1)
{
j=n/(n/i);w[++m]=n/i;
//if(w[m]<=lim) id1[w[m]]=m; else id2[n/w[m]]=m;
h[m]=(w[m]-1)%mod;g[m]=((w[m]+2)%mod)*((w[m]-1)%mod)%mod;
if(g[m]&1)g[m]+=mod;g[m]/=2;
}
for(int j=1;j<=pnum;++j) for(int i=1;i<=m && (ll)pri[j]*pri[j]<=w[i];++i)
{
int k=id(w[i]/pri[j]);//printf("%lld %d %d\n",w[j],i,k);
//cerr<<w[i]<<" "<<pri[j]<<" "<<w[i]/pri[j]<<" "<<k<<endl;
g[i]=(g[i]-(g[k]-sq[j-1])*pri[j]%mod+mod)%mod;
h[i]=((h[i]-h[k]+j-1)%mod+mod)%mod;
}
printf("%lld\n",S(n,1)+1);
return 0;
}