简单介绍

由Min_25在他的博客中介绍的做法，又称Min_25筛。
用于求积性函数 $f(n)$ 的前缀和，其中要求 $f(p)$ 可以表示成多项式，并且 $f(p^k)$ 可以快速算出。
能够在 $O(frac{n^{0.75}}{log})$ 的时间内算出。

基本思路

使用 $DP$ 先求出所有素数的 $f(p)$ 的前缀和，然后再通过枚举最小的质因子计算出合数的前缀和。

求素数的前缀和

下面的所有除法都为整除。
我们之后要求对所有的 $m=frac{n}{i}$ 求 $sum_{isubseteq prime,i<=m}f(i)$
为了方便计算我们将 $f(i)$ 拆成几个完全积性函数的和。使得这些积性函数在素数 $p$ 的位置的和与 $f(i)$ 一致。
接下来考虑一个完全积性函数 $f(i)$ 在素数位置的和怎么求。
那么用容斥，考虑用所有的和减去合数的和。
$g(n,j)=sum_{i=1}^n[isubseteq prime|minp(i)>pri[j]]f(i)$
从 $g(n,j)$ 求 $g(n,j+1)$ 的过程相当于是一个埃氏筛法（即找出 $sqrt n$ 下的质数，将它的倍数筛掉，相当于现在筛到的数的最小因子为当前质数）。
那么如果 $pri[j]^2>n$ ，则 $pri[j]$ 不能继续筛了，所以：
$g(n,j)=g(n,j-1)$ 。
否则 $pri[j]$ 可以筛，即为：
$g(n,j)=g(n,j-1)-pri[j](g(n/pri[j],j-1)-sum_{k<j}f(pri[k]))$
因为 $g$ 里面有小于 $pri[j]$ 的质数，但是它们乘上 $pri[j]$ 的最小因子并不是 $pri[j]$ ，不满足当前筛的数的范围——最小质因子为 $pri[j]$ ，所以减去。
最后求得 $g(n,|P|)$ 即为素数的和了。

亿点点细节

一般使用 $p^k$ 作为完全积性函数，初值 $g(n,0)$ 就是一个自然数幂求和。
这个东西用滚动数组求。
由于要求所有 $g(frac{n}{i},|P|)$ ，查找的时候可以设一个 $id1,id2$ ，如果 $i<=sqrt(n)$ 就直接存在该位置，否则存在 $n/i$ 的位置，可以证明一一对应。然后离散化去储存和转移。
还需要在线筛时预先求出 $sum_{k<j}f(pri[k])$ 。
至此我们就能求出 $sum_{i<=m}f(i)$ ，这里的 $f(i)$ 是题意中的函数。
所有计算都不计算1的贡献，因为筛的时候比较方便

求所有数的前缀和

设 $S(n,j)=sum_{i=1}^n[minp(i)>=pri[j]]f(i)$
$S(n,j)=g(n,|P|)-sum_{k<j}f(pri[k])+sum_{k=j}^{pri[k]^2<=n}sum_{q=1}^{pri[k]^{q+1}<=n}S(n/pri[k]^q,k+1)+f(pri[k]^{q+1})$
前面的质数的和，后面合数的和也是枚举最小的质因子的次幂去筛它，注意最后一项算不到要补上。
出口即 $n<pri[j],n=1$ 返回 $0$ 。
因为整除的数是有限且不重复的，所以 $S(n,j)$ 并不需要像杜教筛一样记忆化。
听说时间是 $O(frac{n^{0.75}}{log})$ 的，最快能过 $10^{12}$ 。
相较杜教筛普适性更高了，但是速度也更慢了。
luoguP3235

#include<cstdio>
#include<cmath>
#include<cstring>
#include<algorithm>
#define maxn 500005
#define ll long long 
#define mo 1000000007
using namespace std;

ll n,sqr,i,j,k,m,inv6;
ll id1[maxn],id2[maxn],w[maxn];
ll bz[maxn],tot,pri[maxn],sg1[maxn],sg2[maxn];
ll s[maxn],g1[maxn],g2[maxn];

ll ksm(ll x,ll y){
	ll s=1;
	for(;y;y/=2,x=x*x%mo) if (y&1)
		s=s*x%mo;
	return s;
}

void getpri(){
	for(i=2;i<=sqr;i++){
		if (!bz[i]) {
			pri[++tot]=i;
			sg1[tot]=(sg1[tot-1]+i)%mo;
			sg2[tot]=(sg2[tot-1]+1ll*i*i)%mo;
		}
		for(j=1;j<=tot&&i*pri[j]<=sqr;j++){
			bz[i*pri[j]]=1;
			if (i%pri[j]==0) break;
		}
	}
}

void prepare(){
	getpri();
	for(i=1;i<=n;i=j+1){
		j=n/(n/i),w[++m]=n/i;
		if (w[m]<=sqr) id1[w[m]]=m; else
			id2[n/w[m]]=m;
		ll tmp=w[m]%mo; 
		g1[m]=(tmp+1)*tmp/2%mo-1;
		g2[m]=tmp*(tmp+1)%mo*(tmp*2+1)%mo*inv6%mo-1;
	}
	for(j=1;j<=tot;j++){
		for(i=1;i<=m&&pri[j]*pri[j]<=w[i];i++){
			int k=(w[i]/pri[j]<=sqr)?id1[w[i]/pri[j]]:id2[n/(w[i]/pri[j])];	
			(g1[i]-=pri[j]*(g1[k]-sg1[j-1]))%=mo;
			(g2[i]-=pri[j]*pri[j]%mo*(g2[k]-sg2[j-1]))%=mo;
		}
	}
}

ll S(ll x,int j){
	if (x==1||pri[j]>x) return 0;
	int k=(x<=sqr)?id1[x]:id2[n/x];
	ll sum=g2[k]-g1[k]-sg2[j-1]+sg1[j-1];
	for(k=j;k<=tot&&pri[k]*pri[k]<=x;k++){
		ll p=pri[k],q=p%mo;
		while (p*pri[k]<=x){
			sum+=S(x/p,k+1)*(q*(q-1)%mo)%mo;
			p=p*pri[k],q=q*pri[k]%mo,sum+=q*(q-1)%mo;
		}
	}
	return sum%mo;
}

int main(){
	scanf("%lld",&n),sqr=sqrt(n),inv6=ksm(6,mo-2);
	prepare();
	printf("%lld",S(n,1)+1);
}