Min_25 筛

用来求积性函数 (f(x)) 的前缀和。要求其在质数 (p) 处的取值为多项式，并且 (f(p^k)) 可以快速计算。

因为多项式可以拆成单项式，所以只用求出形如 (f(p)=p^k) 的前缀和，然后加起来就行。

设 (g(n,j)) 为小于等于 (n) 的所有质数和最小质因子大于第 (j) 个质数的合数的 (k) 次方和，即：

[large g(n,j)=sum_{i=1}^n [i in Prime or ext{minp}(i) > p_j] i^k ]

考虑递推，得：

[largeegin{aligned} g(n,j)&=g(n,j-1)-sum_{i=1}^n[i otin Prime and ext{minp}(i)=p_j]i^k \ &=g(n,j-1)-p_j^ksum_{i=1}^{frac{n}{p_j}}[ ext{minp}(i) geqslant p_j]i^k \ &=g(n,j-1)-p_j^kleft( g(frac{n}{p_j},j-1)-g(p_{j-1},j-1) ight)[p_j^2leqslant n] end{aligned} ]

(g(n,j-1)) 比 (g(n,j)) 多的部分为小于等于 (n) 的最小质因子为 (p_j) 的合数的 (k) 次方和，提出 (p_j) 后，剩下的部分要求最小质因子大于等于 (p_j)，即 (g(frac{n}{p_j},j-1))，但这里多算了质数，所以还要减去 (g(p_{j-1},j-1))。考虑到只有 (p_j leqslant sqrt n) 时，后一项才有值，因此设 (sp(n)=sum_{i=1}^n p_i^k)，得 (g(p_{j-1},j-1)=sp(j-1))，这一部分可以线性筛预处理。

一开始处理出 (g(n,0))，然后递推计算到 (g(n,x))，其中 (x= ext{maxp}(n))，也就是 (g(n,x)) 为小于等于 (n) 的所有质数的 (k) 次方和。

因为 (n) 过大，无法求出所有 (g(n,x))，但是发现形如 (leftlfloor frac{n}{a} ight floor) 的取值只有 (O(sqrt n)) 个，因此只用求这 (O(sqrt n)) 个位置的 (g(n,x)) 即可。

处理出 (g(n,x)) 的复杂度为 (O(frac{n^{frac{3}{4}}}{log n}))。

求前缀和时，类似的设 (s(n,j)) 为小于等于 (n) 的所有最小质因子大于 (p_j) 的数处的函数值的和，即：

[large s(n,j)=sum_{i=1}^n [ ext{minp}(i)>p_j]f(i) ]

将其分为两部分来考虑，分为大于 (p_j) 的质数和最小质因子大于 (p_j) 的合数，质数部分即为 (g(n,x)-sp(j))，合数部分枚举最小质因子即可，得：

[large s(n,j)=g(n,x)-sp(j)+sum_{k>j and p_k^eleqslant n}f(p_k^e)left( s(frac{n}{p_k^e},k)+[e e 1] ight) ]

因为当 (e e 1) 时没有计算 (f(p_k^e)) 的值，所以还要加上 ([e e 1])。实现时就直接递归计算，不用记忆化。

递归的复杂度为 (O(n^{1-epsilon }))。

求积性函数 (f(x)) 的前缀和，其满足 (f(p^k)=p^k(p^k-1))：

#include<bits/stdc++.h>
#define maxn 200010
#define p 1000000007
#define inv2 500000004
#define inv6 166666668
using namespace std;
typedef long long ll;
template<typename T> inline void read(T &x)
{
    x=0;char c=getchar();bool flag=false;
    while(!isdigit(c)){if(c=='-')flag=true;c=getchar();}
    while(isdigit(c)){x=(x<<1)+(x<<3)+(c^48);c=getchar();}
    if(flag)x=-x;
}
ll n,t,tot,cnt;
ll pri[maxn],g1[maxn],g2[maxn],sp1[maxn],sp2[maxn],id1[maxn],id2[maxn],q[maxn];
bool tag[maxn];
int id(ll x)
{
    return x<=t?id1[x]:id2[n/x];
}
void init()
{
    for(int i=2;i<=t;++i)
    {
        if(!tag[i]) pri[++tot]=i;
        for(int j=1;j<=tot;++j)
        {
            int k=i*pri[j];
            if(k>t) break;
            tag[k]=true;
            if(i%pri[j]==0) break;
        }
    }
}
ll s(ll n,int j)
{
    if(pri[j]>=n) return 0;
    int x=id(n);
    ll v=((g2[x]-sp2[j]-(g1[x]-sp1[j]))+p)%p;
    for(int k=j+1;k<=tot&&pri[k]*pri[k]<=n;++k)
    {
        ll now=pri[k];
        for(int e=1;now<=n;++e,now*=pri[k])
        {
            ll val=now%p;
            v=(v+val*(val-1+p)%p*(s(n/now,k)+(e!=1)))%p;
        }
    }
    return v;
}
int main()
{
    read(n),t=sqrt(n),init();
    for(int i=1;i<=tot;++i)
        sp1[i]=(sp1[i-1]+pri[i])%p,sp2[i]=(sp2[i-1]+pri[i]*pri[i])%p;
    for(ll l=1,r;l<=n;l=r+1)
    {
        ll v=n/l;
        r=n/v,q[++cnt]=v,v%=p;
        g1[cnt]=((v+1)*v%p*inv2-1+p)%p;
        g2[cnt]=(v*(v+1)%p*(2*v+1)%p*inv6-1+p)%p;
        if(q[cnt]<=t) id1[q[cnt]]=cnt;
        else id2[n/q[cnt]]=cnt;
    }
    for(int j=1;j<=tot;++j)
    {
        for(int i=1;i<=cnt&&pri[j]*pri[j]<=q[i];++i)
        {
            int x=id(q[i]/pri[j]);
            g1[i]=(g1[i]-pri[j]*(g1[x]-sp1[j-1]+p)+p)%p;
            g2[i]=(g2[i]-pri[j]*pri[j]%p*(g2[x]-sp2[j-1]+p)+p)%p;
        }
    }
    printf("%lld",(s(n,0)+1)%p);
    return 0;
}