Min_25筛学习笔记

太仙了(

Min_25筛是一种快速处理积性函数前缀和问题的算法

使用条件是在$x$是素数时，$f(x)$的表达形式是一个关于$x$的多项式，并且可以快速的计算$f(x^k)$

主要分成三个部分：

第一个部分：素数的函数和

第二个部分：合数的函数和

第三个部分：$1$

先来解决第一部分：

定义函数$g(n,j)$,$g(n,j)$表示的是在$n$以内，素数或者最小质因子大于$Prime_j$的数的函数值的和

注意：此时我们假装所有的函数值都是按照素数来算的。

考虑怎么转移这个东西

首先，从$g(n,j-1)$到$g(n,j)$，少的是最小质因子是$prime_j$的那些合数的贡献。

于是分类讨论一下

如果${prime_j}^{2}>n$，就说明此时最小质因子是$prime_j$的合数已经大于$n$了

于是此时$g(n,j)=g(n,j-1)$

否则的话，因为是积性函数，我们就先把$prime_j$提出来

于是这一部分就是$prime_j * g(\frac{n}{prime_j},j-1)$

但是这里多了点东西

多了什么？

多了素数的函数值和

于是把素数的函数值和加回去即可

这一部分显然是$g(prime_j -1,j-1)$

于是式子就是$g(n,j-1)-prime_j*[g(\frac{n}{prime_j},j-1) - g(prime_j -1,j-1)]$

发现在上面的过程中，其实所有涉及到的值只有根号个，我们筛出这一部分值就好了。

$g$数组其实类似于埃氏筛法第$j$个质数的筛选过程

如果所有数字的求法都和素数一样，那么$g$数组其实就是答案，$g(n,0)$

但显然不是

于是就有了下面的处理：

我们重新定义$F(n,j)$，表示$n$内最小质因子大于等于$prime_j$的和，但是这里这个$F$表示的是积性函数的和

类似的处理

考虑计算$F(n,j)$

素数的部分:$g(n,|P|)-\sum_{1}^{j-1}f(Prime_j)$

合数的部分类似的枚举，这一次我们枚举最小质因数和它的次数

即$\sum_{k>=j} \sum_e f({prime_k}^{e}) * [F(\frac{n}{{prime_k}^{e}},k+1) + f({prime_k}^{e+1})] $

这里需要注意的是，因为${prime_k}^{e}$被筛了，所以要加回去

所以整个式子就是

$g(n,|P|)-\sum_{1}^{j-1}f(Prime_j)$+$\sum_{k>=j} \sum_e f({prime_k}^{e}) * [F(\frac{n}{{prime_k}^{e}},k+1) + f({prime_k}^{e+1})] $

根据这个计算即可。

例题：

1、

LOJ6053 简单的函数

https://loj.ac/p/6053

$f(p)=p-1$

$f(p^k)=p xor k$

维护一个一次的和跟一个$1$的积性函数即可。

#include <bits/stdc++.h>
using namespace std;
long long N;
const long long fish=1e9+7;
long long mx;
const int xx=sqrt(1e11)+100;
long long g[xx],sp[xx],h[xx],Prime[xx],w[xx];
int m,id1[xx],id2[xx],Index;
bool IsPrime[xx];
void Pre(){
    IsPrime[1]=true;
    for (long long i=2;i<=mx;i++){
        if (!IsPrime[i]) Prime[++Index]=i,sp[Index]=(sp[Index-1]+i)%fish;
        for (int j=1;j<=Index&&1ll*i*Prime[j]<=mx;j++){
            IsPrime[1ll*i*Prime[j]]=true;
            if (i%Prime[j]==0) break;
        }
    }
}
long long Pow(long long x,int y){
    long long ans=1;
    while (y){
        if (y&1) ans=1ll*ans*x%fish;
        x=1ll*x*x%fish;
        y>>=1;
    }
    return ans;
}
long long GetAns(long long x,int y){
    if (x<=1||Prime[y]>x) return 0;
    int k=(x<=mx)?id1[x]:id2[N/x];
    long long ans=(g[k]-h[k]-sp[y-1]+y-1+fish)%fish;
    if (y==1) ans+=2;
    for (int i=y;i<=Index&&1ll*Prime[i]*Prime[i]<=x;i++){
        long long t1=Prime[i],t2=1ll*Prime[i]*Prime[i];
        for (long long e=1;t2<=x;e++,t1=t2,t2*=Prime[i]){
            ans=(ans+1ll*(Prime[i]^e)*GetAns(x/t1,i+1)%fish+(Prime[i]^(e+1ll))%fish)%fish;
        }
    }
    return ans;
}
int main(){
    cin>>N;
    mx=sqrt(N);
    Pre();
    long long inv=Pow(2,fish-2);
    for (long long i=1,j;i<=N;i=j+1ll){
        j=N/(N/i);w[++m]=N/i;
        g[m]=1ll*w[m]%fish*1ll*((w[m]+1)%fish)%fish;
        g[m]=1ll*g[m]*inv%fish;
        g[m]=(g[m]-1ll+fish)%fish;
         h[m]=(w[m]-1+fish)%fish;
        if (w[m]<=mx) id1[w[m]]=m;
        else id2[j]=m;
    }
    for (int j=1;j<=Index;j++)
        for (int i=1;i<=m&&1ll*Prime[j]*Prime[j]<=w[i];i++){
            long long k=w[i]/Prime[j];
            if (k<=mx) k=id1[k];else k=id2[N/k];
            (g[i]-=1ll*Prime[j]*(g[k]-sp[j-1])%fish)%=fish;
            (h[i]-=h[k]-j+1)%=fish;
        }    
        //cout<<"qwq"<<endl;
    int ans=GetAns(N,1)+1;
    cout<<(ans+fish)%fish;
    return 0;
}

 

I.Distance

 

下文涉及除法不特殊说明均为下取整

考虑$x={p_1}^{k_1}*{p_2}^{k_2}...{p_n}^{k_n}$

那么显然，$x$能走到$y$的话，一定是沿着质因子一直走

也就是其实我们只需要对质因子分开考虑

枚举每个质因子的贡献

其实是$p_i$*$(i-j)$*(cnt_i - cnt_j)

其中$cnt_i$表示$p_i$指数为$i$的数字的数量

然后进一步考虑优化

考虑$cnt_i$怎么算

容斥一下，$\frac{n}{p^i}$是至少$i$个$p$,$\frac{n}{p^(i+1)}$是至少$i+1$个p

两个对减一下就行了

对于大于$\sqrt(n)$的数字来说，每个数字只会枚举到$0$次和$1$次

也就是$p*cnt_p*cnt_1$

展开一下，即$p*(n-\frac{n}{p})*(n)$

整除分块一下，发现其实要维护的是素数的前缀和

就相当于min_25筛中，$f(p)=p$的函数前缀

维护一下即可

对于小于$\sqrt {n}$的部分暴力即可

#include <bits/stdc++.h>
using namespace std;
long long N;
const long long fish=1e9+7;
long long mx;
const int xx=6e5+10;
long long g[xx],sp[xx],h[xx],Prime[xx],w[xx];
int m,id1[xx],id2[xx],Index;
bool IsPrime[xx];
void Pre(){
    IsPrime[1]=true;
    for (long long i=2;i<=mx;i++){
        if (!IsPrime[i]) Prime[++Index]=i,sp[Index]=(sp[Index-1]+i)%fish;
        for (int j=1;j<=Index&&1ll*i*Prime[j]<=mx;j++){
            IsPrime[1ll*i*Prime[j]]=true;
            if (i%Prime[j]==0) break;
        }
    }
}
long long Pow(long long x,int y){
    long long ans=1;
    while (y){
        if (y&1) ans=1ll*ans*x%fish;
        x=1ll*x*x%fish;
        y>>=1;
    }
    return ans;
}
long long GetAns(long long x,int y){
    if (x<=1||Prime[y]>x) return 0;
    int k=(x<=mx)?id1[x]:id2[N/x];
    long long ans=(g[k]-h[k]-sp[y-1]+y-1+fish)%fish;
    if (y==1) ans+=2;
    for (int i=y;i<=Index&&1ll*Prime[i]*Prime[i]<=x;i++){
        long long t1=Prime[i],t2=1ll*Prime[i]*Prime[i];
        for (long long e=1;t2<=x;e++,t1=t2,t2*=Prime[i]){
            ans=(ans+1ll*(Prime[i]^e)*GetAns(x/t1,i+1)%fish+(Prime[i]^(e+1ll))%fish)%fish;
        }
    }
    return ans;
}
int Get(long long x){
    if (x>mx) return id2[N/x];
    else return id1[x];
}
int main(){
    cin>>N;
    mx=sqrt(N);
    Pre();
    long long inv=Pow(2,fish-2);
    for (long long i=1,j;i<=N;i=j+1ll){
        j=N/(N/i);w[++m]=N/i;
        g[m]=1ll*w[m]%fish*1ll*((w[m]+1)%fish)%fish;
        g[m]=1ll*g[m]*inv%fish;
        g[m]=(g[m]-1ll+fish)%fish;
        if (w[m]<=mx) id1[w[m]]=m;
        else id2[j]=m;
    }
    for (int j=1;j<=Index;j++)
        for (int i=1;i<=m&&1ll*Prime[j]*Prime[j]<=w[i];i++){
            long long k=w[i]/Prime[j];
            if (k<=mx) k=id1[k];else k=id2[N/k];
            (g[i]-=1ll*Prime[j]*(g[k]-sp[j-1])%fish)%=fish;
    }
    long long ans=0,Pre=g[Get(mx)];
    for (long long i=mx;i;i--){
        int k=Get(N/i);
        ans=(ans+1ll*i*(N-i)%fish*((g[k]-Pre+fish)%fish)%fish)%fish;
        Pre=g[k];
    }
    for (int i=1;i<=Index;i++){
        int res=0;
        for (long long j=1,a=Prime[i];a<=N;j++,a*=Prime[i]){
            for (long long k=0,b=1;k<=j-1;k++,b*=Prime[i]){
                ans+=Prime[i]*(j-k)%fish*((N/b-N/b/Prime[i])%fish)%fish*((N/a-N/a/Prime[i])%fish)%fish;
            }
        }
    }
    cout<<ans*2%fish;
    return 0;
}