min25筛

Posted on Edited on In

①概述：

作用：快速求某些积性函数f(x)的前缀和

对f(x)的要求： 要求为积性函数，且：：能快速求出质数处的前缀和：能快速求出的值

由于要求1，故绝大多数能用min25筛的情况下，f(p)都是关于p的简单多项式（简单多项式的前缀和可以快速求出）

比如： 令为的前缀和，即若则若则

复杂度：

②前置知识：

1.埃氏筛

分析一波埃氏筛复杂度：由素数定理，设表示不大于的数中有多少个质数，则可知：是素数的概率为复杂度我们积分近似一下：设，有所以埃氏筛的复杂度为

2.线性筛

复杂度线性 if(i%p[j]==0)break;

i是p[j]的倍数时，i*p[j+1]这个合数肯定也能被p[j]筛掉，所以不用继续筛下去了

（因为i中含有p[j]，而p[j]肯定小于p[j+1]）

③正式开始：

求

part1:质数部分 （dp）

需要引入 完全积性函数 F(x)，且满足质数处有F(p)=f(p) 令表示第个质数是质数的最小质因子质数部分的和即为因为以内，没有最小质因子的合数 g(n,j)，如何状态转移？

其实就是考虑第j轮埃氏筛，会筛掉哪些数！

g(30,1) : 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30

g(30,2) : 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30

g(30,3) : 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30

g(n,j-1) -> g(n,j)，减小的值为最小质因子为p[j]的合数的值！

设被筛掉的合数为x,有: 用到的是完全积性函数的性质

可得出转移方程当当

注意以上的是我们引入的完全积性函数，仅满足在质数处和相同是递推起点，即在质数处的前缀和，这个需要能够快速求出（即之前所说的要求）第二部分要用到的，是即在质数处的前缀和

part2:加上合数部分的贡献，求出答案

质数处的和已经在part1求出。

容易想到，可以通过枚举最小质因子的幂次的方式，来统计合数部分的贡献（且只用到积性函数的性质就够了） 具体地，令的最小质因子则有 再简单解释一下： 质数部分对的贡献即为合数部分，考虑任何正整数都能表示成我们枚举最小质因子，然后就转化成了子问题。由于积性函数的性质，和直接相乘即可（因为我们保证和是互质的）的意义：这个数就是质数的幂次不能算，因为统计的是合数的贡献

④复杂度分析

part1:(dp模拟埃氏筛)

把放缩到，把放缩到（根号增长远大于对数）

part2:

写法：递归版（无需记忆化），渐进意义下是线性的，但是时跑得很快

（证明见集训队论文《一些特殊的数论函数求和问题》朱震霆）

写法：又名洲阁筛，递推求，复杂度同

⑤代码实现

洛谷模板题：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const ll mod=1e9+7;
ll ksm(ll a,ll b)
{
    ll res=1;
    while(b)
    {
        if(b&1)res=res*a%mod;
        a=a*a%mod , b>>=1;
    }
    return res;
}
ll ny(ll x){return ksm(x,mod-2);}


namespace MIN25//f1(x)=x , f2(x)=x*x
{
    const ll N=1e10+7 , __N=1e5+5;
    ll n,__n;//筛到__n(根号n)即可 
    //因为n以内,没有 最小质因子>根号n 的合数
    
    void read_n(){scanf("%lld",&n) , __n=sqrt(n);}
    
    bool vis[__N];
    ll p[__N];int jsp;//根号n以内的p的个数 
    ll sp1[__N/*jsp*/],sp2[__N/*jsp*/];//F在质数处的前缀和 
    void xxs()//线性筛预处理
    {
        for(int i=2;i<=__n;++i)
        {
            if(!vis[i])
            {
                p[++jsp]=i; 
                sp1[jsp] = (sp1[jsp-1]+i)%mod;
                sp2[jsp] = (sp2[jsp-1]+1LL*i*i)%mod;
            }
            
            for(int j=1;j<=jsp && p[j]*i<=__n;++j)
            {
                vis[p[j]*i]=true;
                if(i%p[j]==0)break;
            }
        }
    }
    
    ll w[__N*2/*tot*/];int tot;//所有整除分块的值（返回n/l），tot为块号 
    int id1[__N*2/*n/l的值*/],id2[__N*2/*相当于(n/l)的值*/];//w的逆映射，返回块号tot 
    inline int get_id(ll val)//val是n的整除分块的值（n/l） 
    {
        if(val<=__n)return id1[val];
        else return id2[n/val];
    }
    ll g1[__N*2/*tot*/],g2[__N*2/*tot*/];//初始为g(n,0) (对n的每个整除分块的值都要算g的)
    
    void work_g()
    {
        for(int j=1;j<=jsp;++j)//g(n/... , j)，只保存n/...这一维
        {
            for(int i=1;i<=tot && p[j]*p[j]<=w[i];++i)//每次要修改所有n的整除分块的g
            //注意顺序，w数组的值是递减的！ 
            {
                ll tmp=w[i]/p[j];
                int tmptot = get_id(tmp);
                
                g1[i] -= p[j]/*f1*/*(g1[tmptot]-sp1[j-1]); 
                g1[i]%=mod;
                g2[i] -= p[j]*p[j]/*f2*/%mod*(g2[tmptot]-sp2[j-1]);
                g2[i]%=mod;
            }
        }
    }
    ll S(ll i,int j)//i是n/... (的值) 
    {
        if(p[j]>=i)return 0;
        int itot = get_id(i);
        ll res = (g2[itot]-g1[itot]) - (sp2[j]-sp1[j]);//质数部分
        res%=mod;
        
        for(int k=j+1;p[k]*p[k]<=i && k<=jsp;++k)
        {
            ll pe=p[k];//p[k]的e次方
            for(int e=1;pe<=i;++e,pe=pe*p[k])//pe这里先别取模！要作为循环条件！
            {
                ll x=pe%mod;
                res += (x*(x-1)%mod)/*f(pe)*/ * (S(i/pe,k)+(e>1));
                res%=mod;
            }
        }
        
        return res;
    }
    
    ll work__min25() 
    {
        read_n() , xxs();
        
        ll inv6=ny(6);
        auto calc = [inv6](ll i,int op){
            i%=mod;
            if(op==1)return i*(i+1)/2%mod;
            return i*(i+1)%mod*(i+i+1)%mod*inv6%mod;
        };
        
        for(ll l=1,r;l<=n;l=r+1)
        {
            r=min(n,n/(n/l));
            
            w[++tot]=n/l;
            if(w[tot]<=__n)id1[w[tot]] = tot;
            else id2[n/w[tot]] = tot;
            
            g1[tot]=calc(w[tot],1)-1 , g2[tot]=calc(w[tot],2)-1;
            //g1,g2记得减掉F(1)的值 ！ 
        }
        
        work_g();
        ll res=S(n,0)+1;
        res=(res%mod+mod)%mod;
        return res;
    }
}

int main()
{
    printf("%lld\n",MIN25::work__min25());
    return 0;
}