Min_25筛
参考
题意：给定一积性函数 \(f\)，其在 \(p^c\) 处的取值是一个关于 \(p^c\) 的低次多项式，求 \(\sum_{i=1}^n f(i)\)。
\(1\le n\le 10^{10}\)
先处理质数处的取值。
拆质数处的多项式，依次处理每个单项式。
考虑把和分为两部分，质数和非质数处的取值。
考虑一积性函数 \(g\)，其在质数处的取值和 \(f\) 对应的单项式相同。
那么 \(\sum_{i=1}^n g(i)\) 是容易求的，我们考虑去掉其中的合数，可以用每个合数的最小值因子筛掉他。
设 \(h[i,j]\) 表示 \(1...i\) 中筛去了最小质因子 \(\in\) 前 \(j\) 个质数的合数后，剩下数取值的和。
记 \(p_i\) 为第 \(i\) 个质数。考虑从 \(h[i,j-1]\) 转移，我们需要去掉最小质因子为 \(p_j\) 的数。这些数除以 \(p_j\) 后，其最小值因子都 \(\ge p_j\)。
不难发现，最小值因子为 \(p_j\) 的数的取值总和为 \(g(p_j)\cdot(h[\lfloor \frac n{p_j}\rfloor,j-1]-\sum_{k=1}^{j-1}g(k))\)
那么
注意一点，就是 \(h[i,j]\) 始终保留 \(g(1)\) 的值，因为 \(1\) 永远不会被筛掉。
再考虑合数，下文默认 \(h[i,j]\) 为已考虑的所有单项式的和。
一种想法是通过枚举他的最小值因子来求。
设 \(S[i,j]\) 表示 \(2...i\) 中筛去最小值因子 \(\in\) 前 \(j\) 个质数后剩下的数的 \(f\) 取值之和。
首先算上质数，贡献为 \(h[i,m]\)（\(m\) 为 \(\sqrt n\) 内质数个数）。枚举最小值因子为第 \(k\) 个质数，其次数为 \(c\)，贡献为 \(f(p_k^c)\cdot (S[\lfloor \frac n{p_k^c} \rfloor,k]+[c>1])\)，其中 \(c>1\) 是因为此时 \(p_k^c\) 自身也是一个合数。
Q: 为什么质数要单独分开？
A: 因为质数值域很大，而合数的最小值因子 \(\le 10^5\)
时间复杂度不会证，\(O(\frac {n^{\frac 34}}{\log n})\)。
答案是 \(S[n,0]-\sum g(1)+1\)，注意 \(S[n,0]\) 算上的是 \(1\) 在 \(g\) 处的单项式总和，而不是 \(f(1)\)。
#include<bits/stdc++.h>
#define pr pair<ll,ll>
#define x first
#define y second
#define mkp(a,b) make_pair(a,b)
#define pb push_back
#define ll long long
using namespace std;
const ll maxn=2e6+10,mod=1e9+7;
ll n,sq,pri[maxn],m,id1[maxn],id2[maxn],sum1[maxn],sum2[maxn],w[maxn],len,g1[maxn],g2[maxn];
bool b[maxn];
void xxs()
{
sum1[0]=sum2[0]=1;
for(ll i=2;i<=sq;i++)
{
if(!b[i])
{
pri[++m]=i;
sum1[m]=sum1[m-1]; sum2[m]=sum2[m-1];
sum1[m]+=i, sum2[m]+=i*i;
sum1[m]%=mod; sum2[m]%=mod;
}
for(ll j=1;j<=m&&i*pri[j]<=sq;j++)
{
ll k=i*pri[j]; b[k]=true;
if(i%pri[j]==0) break;
}
}
}
ll id(ll x)
{
if(x<=sq) return id1[x];
return id2[n/x];
}
ll Sum1(ll n)
{
n%=mod;
return n*(n+1)%mod*(mod/2+1)%mod;
}
ll Sum2(ll n)
{
n%=mod;
return n*(n+1)%mod*(2*n+1)%mod*166666668%mod;
}
ll F(ll n,ll k)
{
if(pri[k]>n) return 0;
ll ans=((g2[id(n)]-g1[id(n)]+mod)%mod-(sum2[k]-sum1[k]+mod)%mod+mod)%mod;
for(ll i=k+1;i<=m&&pri[i]*pri[i]<=n;i++)
for(ll x=1,y=pri[i];y<=n;y*=pri[i],++x)
{
ans=(ans+y%mod*(y%mod+mod-1)%mod*(F(n/y,i)+(x>1)))%mod;
}
return ans;
}
int main()
{
scanf("%lld",&n); sq=sqrt(n);
xxs();
for(ll i=1;i<=n;i++)
{
ll d=n/i, r=n/d; w[++len]=d;
g1[len]=Sum1(d); g2[len]=Sum2(d);
if(d<=sq) id1[d]=len;
else id2[n/d]=len;
i=r;
}
for(ll i=1;i<=m;i++)
{
for(ll j=1;j<=len&&pri[i]*pri[i]<=w[j];j++)
{
g1[j]=(g1[j]-pri[i]*(g1[id(w[j]/pri[i])]-sum1[i-1])%mod+mod)%mod;
g2[j]=(g2[j]-pri[i]*pri[i]%mod*(g2[id(w[j]/pri[i])]-sum2[i-1])%mod+mod)%mod;
}
}
printf("%lld",(F(n,0)+1)%mod);
return 0;
}
