MIN_25 筛
适用范围
亚线性 (\(O(\dfrac{n^\frac{3}{4}}{\log n})\)) 求解部分积性函数的前缀和。
要求 \(f(x)\) 是关于 \(x\) 的多项式函数且项数不是很多,次数不是很大，而且可以较快地求出 \(f(p^c)\),其中 \(p\) 是一个素数
描述
给定 \(n\) 和积性函数 \(f\)
求 \(\sum\limits_{i=1}^{n} f(i)\)
符号规定：
-
默认 \(p\in \mathbb{P}\) 其中 \(\mathbb{P}\) 为素数集。
-
\(minp(x)\) 表示 \(x\) 的最小质因子。
-
\(p_i\) 表示第 \(i\) 小的质数，特别地，\(p_0=1\)
-
\(S_k(x)=\sum\limits_{i=2}^x [p_k\leq minp(i)] f(i)\)
显然 \(\sum\limits_{i=1}^{n} f(i)=S_1(n)+f(1)\)
考虑如何求 \(S_k(x)=\sum\limits_{i=2}^x [p_k\leq minp(i)] f(i)\)
枚举每个 \(p_j,k\leq j\)
\(S_k(x)=\sum\limits_{k\leq j,p_j\leq x} \sum\limits_{1\leq c,{p_j}^c\leq x} f({p_j}^c)(f(1)+S_{j+1}(\lfloor\dfrac{x}{{p_j}^c}\rfloor))\)
有大量的 \(p_j\) 的指数 \(c\) 只能取到 \(1\) 且 \(S_{j+1}(\lfloor\dfrac{x}{{p_j}^c}\rfloor)=0\) 于是把和式改写成如下形式
\(S_k(x)=\sum\limits_{k\leq j,{p_j}^2\leq x}\sum\limits_{1\leq c,{p_j}^c\leq x}f({p_j}^c)(f(1)[c>1]+S_{j+1}(\lfloor\dfrac{x}{{p_j}^c}\rfloor)) + f(1) \sum\limits_{k\leq j,p_j\leq x} f(p)\)
记 \(F(x)=\sum\limits_{p\leq x}f(p)\)
则 \(S_k(x)=\sum\limits_{k\leq j,{p_j}^2\leq x}\sum\limits_{1\leq c,{p_j}^c\leq x}f({p_j}^c)(f(1)[c>1]+S_{j+1}(\lfloor\dfrac{x}{{p_j}^c}\rfloor)) + f(1)(F(x)-F({p_k}-1))\)
如果我们能求出 \(F\) 在上式中的取值，可以证明计算 \(S_1(n)\) 的复杂度是 \(O(\dfrac{n^\frac{3}{4}}{\log n})\)
下面来求 \(F\) 在上式的取值。
/*
由于 \(x\) 是合数时 \(minp(x)\leq \sqrt{x}\),当 \(k>\sqrt{n}\) 时，\(S_k(n)=0\)，所以 \(k\) 的取值有 \(O(\sqrt{n})\) 种
又由数论分块可知,求解 \(S_1(n)\) 的过程中，要用到的 \(S_k(x)\) 的不同的 \(x\) 的值有 \(O(\sqrt{n})\) 种
*/
由于 \(F(x)=\sum\limits_{p\leq x}f(p)\) 且 \(f(p)\) 是关于 \(p\) 的项数较少的多项式，所以可以把每一次项拆开计算。
即 \(F(x)=\sum a_i \sum\limits_{p\leq x}p^{c_i}\)
于是问题变成了求 \(\sum\limits_{p\leq x} p^m\)
要计算 \([2,x]\) 内所有素数的贡献，考虑一个类似埃氏筛的过程。
定义函数 \(g(x,k)=\sum\limits_{i=2}^x [i\in \mathbb{P}~ \texttt{or} ~ p_k<minp(i) ]~i^m\)
也就是埃氏筛筛完前 \(k\) 轮后剩下的数的函数值
由于对于合数 \(x\)， \(minp(x)\leq \sqrt{x}\),所以 \(\sum\limits_{p\leq x}p^m=g(x,\lfloor\sqrt{x}\rfloor)\)
因为 \(m\) 不是很大，所以我们有办法求出 \(g(x,0)=\sum\limits_{i=2}^x i^m\)
然后我们通过 dp 求其他要用的 \(g\) 的值。
考虑 \(g(x,k)\) 从 \(g(x,k-1)\) 转移过来。
我们要去掉第 \(k-1\) 轮完成后剩下的最小质因数是 \(p_k\) 的数。
于是要减去 \({p_k}^mg(\lfloor\dfrac{x}{p_k}\rfloor,k-1)\)
但是对于形如 \(y={p_k}\times q,q\in \mathbb{P} ,q<p\) 的数，我们多减了
于是还需加上 \([1,p_{k-1}]\) 内素数带来的转移的贡献。
综上，\(g(x,k)=g(x,k-1)-p_{k}^mg(\lfloor\dfrac{x}{p_k}\rfloor,k-1)+p_{k}^mg(p_{k-1},k-1)\)
#include<bits/stdc++.h>
using namespace std;
const int mod=1e9+7,N=1e6,inv2=500000004,inv6=166666668;
typedef long long ll;
ll val[N],cnt,n,T,g1[2][N],g2[2][N];
unordered_map<ll,int> pos;
vector<ll> primes;
bitset<N> st;
inline ll getF(ll u){
return (g2[0][pos[u]]-g1[0][pos[u]]+mod)%mod;
}
void sieve(){
ll T=sqrt(n)*2;
for(int i=2; i<=T; i++){
if(!st[i]){
primes.emplace_back(i);
}
for(auto j:primes){
if(i*j>T) break;
st[i*j]=1;
if(i%j==0) break;
}
}
primes.insert(primes.begin(),1);
for(int i=1; i<=cnt; i++){
g1[0][i]=val[i]%mod*(val[i]%mod+1)%mod*inv2%mod;
g2[0][i]=val[i]%mod*(val[i]%mod+1)%mod*(2ll*val[i]%mod+1)%mod*inv6%mod;
}
int s=cnt,x=1;
for(int k=1; k<primes.size(); k++){
x=k&1;
for(int i=s; i; i--){
// cout<<val[i]<<' ';
g1[x][i]=g1[x^1][i];
g2[x][i]=g2[x^1][i];
if(primes[k]*primes[k]>val[i]){
s=i-1;
continue;
}
g1[x][i]=(g1[x][i]-primes[k]*g1[x^1][pos[val[i]/primes[k]]]%mod);
g1[x][i]=(g1[x][i]+primes[k]*g1[x^1][pos[primes[k-1]]]%mod)%mod;
g2[x][i]=(g2[x][i]-primes[k]*primes[k]%mod*g2[x^1][pos[val[i]/primes[k]]]%mod);
g2[x][i]=(g2[x][i]+primes[k]*primes[k]%mod*g2[x^1][pos[primes[k-1]]]%mod)%mod;
}
// cout<<"using prime number "<<primes[k]<<endl;
}
for(int i=s; i; i--){
g1[x^1][i]=g1[x][i];
g2[x^1][i]=g2[x][i];
}
}
ll F(ll k,ll u){
// cout<<k<<' '<<val[u]<<endl;
if(k>=primes.size()) return 0;
if(primes[k]>val[u]) return 0;
ll sum=0;
for(ll j=k; j<primes.size()&&primes[j]*primes[j]<=val[u]; j++){
ll t=1;
for(ll c=1; ;c++){
t=t*primes[j];
if(t>val[u]) break;
sum=(sum+t%mod*(t%mod-1)%mod*((c>1?1:0)+F(j+1,pos[val[u]/t]))%mod);
}
sum%=mod;
}
return ((sum+getF(val[u])-getF(primes[k]-1))%mod+mod)%mod;
}
int main(){
cin>>n;
ll l,r;
for(l=1;l<=n;l=r+1){
r=n/(n/l);
val[++cnt]=n/l;
pos[n/l]=cnt;
// cout<<n/l<<' ';
}
// cout<<endl;
sieve();
cout<<(F(1,1)+1)%mod<<endl;
return 0;
}