Min_25筛
简介
min_25筛据说可以在\(o(\frac{n^{\frac{3}{4}}}{logn})\)处理出含有以下性质的函数f的前缀和:
1.\(f(ab)=f(a)f(b)\),即f是一个积性函数
2.\(f(p^k)\)可以快速计算。
ps:本文没有关于复杂度的证明。。。
预处理
首先要预处理两个东西,一个是\(\sqrt{n}\)(n为询问的值域)内的质数。直接线性筛就好了。用\(pri[i]\)表示第i个质数。设共有\(m\)个质数
另一个是\(g(n,j)\),表示所有\(x\in[1,n]\)中满足x最小质因子大于\(pri[j]\)或者x是质数的\(f(x)\)之和。
这样一来,\(g(n,m)\)表示的就是\([1,n]\)中所有质数的f值之和。这个东西后面会用到。
那么来看一下这个\(g\)值该如何求。
显然,如果\(pri_j^2>n\)那么\(g(n,j)=g(n,j-1)\)。因为这时的\(g(n,j-1)\)已经只表示\([1,n]\)中所有质数的f之和。\(g(n,j)\)并不会比\(g(n,j-1)\)多删除掉任何东西。
如果\(pri_j^2\le n\)呢?我们可以理解为埃氏筛法的过程,\(g(n,j)\)与\(g(n,j-1)\)的差别就是筛掉了\(pri_j\)的倍数。那么就好像可以转移了。问题就在于如何计算出所有\(pri_j\)的倍数所产生的贡献。前面说到这是一个积性函数,所以我们将这些要删除的数全都提出来一个\(pri_j\),那么剩下的就是\([1,\frac{n}{pri_j}]\)了。因为需要\(pri_j\)是这些数字的最小质因子,所以实际上区间\([1,pri_j-1]\)内的数字是不可以的,所以要删除的区间实际上是\([pri_j,\frac{n}{pri_j}]\)所以要删除的数字就是\(f(pri_j)[g(\frac{n}{pri_j},j-1)-g(pri_j-1,j-1)]\)。也就是说\(g(n,j)=g(n,j-1)-f(pri_j)[g(\frac{n}{pri_j},j-1)-g(pri_j-1,j-1)]\)
因为最终我们需要的只有\(g(\lfloor\frac{n}{i}\rfloor,m),1\in [1,n]\)。所以空间只需要开一维,每次处理复杂度是\(o(m)\)的(实际上并不到),类似于整除分块,我们知道\(\frac{n}{i}\)只有\(\sqrt{n}\)级别种取值。复杂度据说是\(o(\frac{n^{\frac{3}{4}}}{logn})\)。
计算答案
上面的东西预处理完了,那么有什么用呢??
我们再定义一个函数\(s(n,j)\)表示\(x\in[1,n]\)中满足\(x\)的最小质因子大于等于\(pri_j\)的\(f(x)\)之和。
最终我们要求的答案就是\(s(n,1)+f(1)\)
上面说到,\(g(n,m)(pri_m^2>n)\)可以表示\([1,n]\)中所有质数的\(f\)值之和。
所以我们将\(s(n,j)\)分为质数和合数两块来处理。
质数的一块显然就是\(g(n,m)-\sum\limits_{k=1}^{j-1}f(pri_k)\)。为什么要减掉后面这一块??因为小于\(pri_j\)的质数不包含在\(s(n,j)\)里面呀~
然后考虑合数的一块该如何求,我们枚举一下这些合数的最小质因子\(k\in[pri_j,pri_m]\)和\(k\)的指数\(e\)。于上方求g的方法类似的,我们可以提出来一个\(pri_k^e\),那么剩下的就是\([1,\frac{n}{pri_k^e}]\),他们的f之和就是\(s(\frac{n}{pri_k^e},k)\)。发现这样无法转移,那么我们只好从\(s(\frac{n}{pri_k^e},k+1)\)转移过来,但是这样\(f(pri_k^{e+1})\)就没被计算,单独加上就好了。
综上所述,
\[s(n,j)=g(n,m)-\sum\limits_{k=1}^{j-1}f(pri_k)+\sum\limits_{k=j}^m\sum\limits_{e=1}^{pri_k^{e+1}\le n}(f(pri_k^{e+1})+f(pri_k^e)s(\frac{n}{pri_k^e},k+1)\]
然后递归计算即可。这里的复杂度据说也是\(o(\frac{n^{\frac{3}{4}}}{logn})\)。
经典例题
定义函数\(f(x)\)满足以下性质,
1.\(f(1)=1\)
2.\(f(p^c)=p\otimes c\)(p为质数)
3.\(f(ab)=f(a)f(b)(a,b互质)\)
求\(\sum\limits_{i=1}^nf(i)\)。\(n\le 10^{10}\)
思路
发现这些性质恰好吻合了我们一开始要求的性质。
发现除2外所有的质数均为奇数,所以就有\(f(p)=p-1\)(p为奇质数)。然后发现这个东西并不是积性函数,没法预处理g了。怎么办?
那就把它拆开,拆成\(f_1(p)=p,f_2(p)=1\),那么就有\(f(p)=f_1(p)-f_2(p)\)。然后按照上述方法分别预处理除关于\(f_1\)的\(g(n,m)\)。关于\(f_2\)的\(h(n,m)\)。
要说明的是,我们一开始将所有的合数全都当成奇质数来处理,因为最后都要“筛”掉的,所以没有影响。
具体细节:
预处理g时,有个式子是\(g(pri_j-1,j-1)\),也就是前\(j-1\)个质数的前缀和。所以预处理质数时同时预处理一下前缀和。
code
/* * @author: wxyww * @date: 2019-12-22 17:42:00 * @last modified time: 2019-12-24 21:58:42 */ #include<cstdio> #include<iostream> #include<cstdlib> #include<cstring> #include<cmath> #include<algorithm> #include<queue> #include<vector> #include<ctime> using namespace std; typedef long long ll; const int n = 2000010,mod = 1e9 + 7; ll read() { ll x = 0,f = 1;char c = getchar(); while(c < '0' || c > '9') { if(c == '-') f = -1; c = getchar(); } while(c >= '0' && c <= '9') { x = x * 10 + c - '0'; c = getchar(); } return x * f; } int tot; ll n,m,w[n],pri[n],sum[n],js,vis[n],g[n],h[n]; void pre() {//筛出质数 for(int i = 2;i <= m;++i) { if(!vis[i]) { pri[++js] = i; sum[js] = sum[js - 1] + i; sum[js] %= mod; } for(int j = 1;j <= js && i * pri[j] <= m;++j) { vis[pri[j] * i] = 1; if(i % pri[j] == 0) break; } } } int id1[n],id2[n]; ll s(ll now,int x) { if(now <= 1 || pri[x] > now) return 0; // printf("%lld %d\n",now,x); int k; if(now <= m) k = id1[now]; else k = id2[n / now]; ll ret = (g[k] - h[k] - sum[x - 1] + x - 1) % mod; if(x == 1) ret += 2;//f(2)当作1计算,实际上f(2)=3 for(int k = x;k <= js && pri[k] * pri[k] <= now;++k) { ll p = pri[k]; for(int e = 1;p * pri[k] <= now;p = p * pri[k],++e) { ret += (pri[k] ^ e) * s(now / p,k + 1) % mod + (pri[k] ^ (e + 1)); ret %= mod; } } return ret; } int main() { // freopen("1.in","w",stdout); n = read(); m = sqrt(n); pre(); // puts("!!!"); for(ll l = 1,r;l <= n;l = r + 1) { r = n / (n / l); ll tmp = n / l; w[++tot] = tmp; if(tmp <= m) id1[tmp] = tot;//数组不够大,通过id1和id2来映射到sqrt(n)以内 else id2[n / tmp] = tot; g[tot] = (tmp + 2) % mod * ((tmp - 1) % mod) % mod; if(g[tot] & 1) g[tot] += mod; g[tot] /= 2; // g[tot] %= mod; h[tot] = tmp - 1; } // for(int i = 1;i <= tot;++i) printf("%lld ",g[i]); for(int j = 1;j <= js;++j) { for(int i = 1;i <= tot && pri[j] * pri[j] <= w[i];++i) {//枚举顺序不能颠倒 ll tmp = w[i] / pri[j]; int k; if(tmp <= m) k = id1[tmp]; else k = id2[n / tmp]; g[i] -= pri[j] * (g[k] - sum[j - 1]) % mod;//枚举顺序不能颠倒的原因 g[i] %= mod; h[i] -= (h[k] - (j - 1)); h[i] %= mod; } } // cout<<g[tot - 2]<<endl; cout<<(s(n,1) + 1 + mod) % mod;//单独把1算上 return 0; }
上一篇: 策划一个爆粉抖音号,你要做的3个布局