Min_25 筛

简介

在开始之前，读者需要了解一些基础的数论内容与关于积性函数的有关定义与基本筛法。

我们定义 \(p_i\) 表示从小到大第 \(i\) 个质数。

不妨设我们要求 \(f\) 函数的前缀和，那么令 \(g(n)\) 表示 \(\sum_{p\le n}f(p)\)，其中 \(p\) 是质数。

然后，不妨假设 \(f\) 是完全积性函数，如果不是就拆成若干项完全积性函数再相加即可。

但是你发现，\(g\) 完全不好求，所以我们考虑去 dp。设 \(g(n,j)=\sum_{i=1}^{n}f(i)[i\in S_j]\)，其中 \(S_j\) 表示质数和最小质因子 \(>p_j\) 的合数的并集。

显然，对于任意合数 \(x\)，一定存在它的一个质因子 \(q\) 使得 \(q\le \sqrt x\)。所以我们要求的就是 \(g(n,k)\)，其中 \(k\) 是最小的数，使得 \(p_k\ge \sqrt n\)。

然后我们考虑从 \(j-1\) 变为 \(j\)。那么显然最小质因子恰好为 \(p_j\) 的合数都会被筛掉。所以有转移：

\[g(n,j)=g(n,j-1)-f(p_j)(g(\lfloor\dfrac{n}{p_j}\rfloor,j-1)-g(p_{j-1},j-1)) \]

因为 \(f\) 是一个完全积性函数，所以我们可以直接把它提出来。然后括号里的前者表示你考虑所有 \(p_j\) 的倍数，去掉一个 \(p_j\) 之后最小质因子 \(>p_{j-1}\) 也就是最小质因子恰好是 \(p_j\) 的所有数的贡献。但是这样直接减掉，会把全体质数减掉，显然之前有一些小质数已经被干掉了，这里再干掉一次就爆了，所以要把这些 \(\le p_{j-1}\) 的小质数加回来。

然后你不太困难的发现，你用到的 \(g(i,j)\) 中的 \(i\) 不会太多，用整除分块之类的状物分析一下可以发现只有大约 \(2\sqrt n\) 个。然后你把这些数预处理出来搞个离散化之类的东西就行了（更准确的说法是重标号）。

显然我们还需要预处理 \(\sqrt n\) 以内的所有质数。

欸但是这里就有问题了，我并不是只用了除一个质数下取整状物啊，我还有小于根号的所有质数啊。

这里细想一下可以发现，前面那 \(2\sqrt n\) 个数就能包括这些数了，证明这里懒得写了。所以是不需要对这些东西做额外操作的。

时间复杂度是 \(O(\dfrac{n^{\frac{3}{4}}}{\log n})\) 的，但是我仍然不知道怎么证明。

然后我们终于做完了 \(g\) 的处理。

现在的目标变成了求 \(s(n)=\sum_{i\le n}f(i)\)。然后，我们做类似 \(g\) 的处理，设 \(s(n,j)=\sum_{i=1}^{n}f(i)[i\in T_j]\)，这里 \(T_j\) 指所有最质因子 \(>p_j\) 的数的集合。

这里的 \(f\) 就不要求是完全积性函数了，只要是积性函数就行，所以我们不用拆了。

然后你考虑，把所有数的贡献拆成质数和合数两部分（当然还有 \(1\)，但是 \(f(1)\) 是 naive 的），所以有：

\[s(n,j)=g(n)-g(p_j)+\sum_{j<k,p_k\le \sqrt n,1\le e,p_k^e\le n}f(p_k^e)(s(\lfloor\dfrac{n}{p_k^e}\rfloor,k)+[e\ne 1]) \]

前半部分显然，就是在算质数的贡献。然后后半部分是在枚举合数的最小质因子以及其次数，然后因为你这个 \(f\) 是积性函数，且括号里面的 \(s\) 里面所有数都与它互质，所以可以直接拆出来。然后这个艾弗森括号是因为如果 \(e=1\)，那么拆出来的是个质数，不应该在这里算；否则应该在这里算。

然后这个东西你直接暴力递归就是对的，甚至连记忆化都不需要。

模板

模板

直接把原函数拆成 \(p^2-p\) 直接上 Min_25 筛即可。

下面是代码：

点击查看代码

#include<bits/stdc++.h>
#define int long long
#define N 200005
#define mod 1000000007
#define pii pair<int,int>
#define x first
#define y second
#define pct __builtin_popcount
#define mpi make_pair
#define inf 2e18
using namespace std;
int T=1,n,id1[N],id2[N],m;
int pre[N],cnt,g1[N],g2[N],v[N];
bool st[N];
int id(int x){
    if(x<N)return id1[x];
    else return id2[n/x];
}
int s1(int x){
    int inv=(mod+1)/2;
    x%=mod;
    return inv*x%mod*(x+1)%mod;
}
int s2(int x){
    int inv=(mod+1)/6;
    x%=mod;
    return inv*x%mod*(x+1)%mod*(x+x+1)%mod;
}
int pw(int x){
    x%=mod;
    return x*x%mod;
}
int f(int x){
    x%=mod;
    return (pw(x)-x+mod)%mod;
}
void init(int n){
    for(int i=2;i<=n;i++){
        if(!st[i])pre[++cnt]=i;
        for(int j=1;j<=cnt&&pre[j]*i<=n;j++){
            st[i*pre[j]]=1;
            if(i%pre[j]==0)break;
        }
    }
}
int s(int x,int y){
    if(pre[y]>=x)return 0;
    int res=(g2[id(x)]-g1[id(x)]-g2[id(pre[y])]+g1[id(pre[y])]+mod*2)%mod;
    for(int i=y+1;i<=cnt&&pre[i]*pre[i]<=x;i++){
        int w=pre[i];
        for(int j=1;w<=x/pre[i];j++,w*=pre[i]){
            (res+=f(w)*s(x/w,i)%mod+f(w*pre[i]))%=mod;
        }
    }
    return res;
}
void solve(int cs){
    if(!cs)return;
    cin>>n;
    init(sqrtl(n)+1);
    for(int l=1,r;l<=n;l=r+1){
        r=n/(n/l);
        v[++m]=n/l;
        if(v[m]<N)id1[v[m]]=m;
        else id2[n/v[m]]=m;
        g1[m]=(s1(v[m])-1+mod)%mod;
        g2[m]=(s2(v[m])-1+mod)%mod;
    }
    for(int j=1;j<=cnt;j++){
        for(int i=1;i<=m&&pre[j]<=v[i]/pre[j];i++){
            g1[i]=(g1[i]-pre[j]*(g1[id(v[i]/pre[j])]-g1[id(pre[j-1])])%mod+mod)%mod;
            g2[i]=(g2[i]-pw(pre[j])*(g2[id(v[i]/pre[j])]-g2[id(pre[j-1])])%mod+mod)%mod;
        }
    }
    cout<<(s(n,0)+1)%mod;
}
signed main(){
	ios::sync_with_stdio(0);
	cin.tie(0);cout.tie(0);
	// cin>>T;
	// init();
	for(int cs=1;cs<=T;cs++){
		solve(cs);
	}
    // cerr<<clock()*1.0/CLOCKS_PER_SEC<<'\n';
	return 0;
}