Min_25筛学习笔记
为方便，下文中的 \(p\) 一律表示质数。
Min_25太神力！
一开始学起来非常费劲比杜教筛难学很多，因为中间的思路比较跳跃，如果没有在一开始理清楚思路的话是很难看明白在干啥的。
Min_25筛的作用
可以在亚线性（有说是 \(O(n^{1-\epsilon})\) 的，也有说是 \(O(\frac{n^{\frac{3}{4}}}{\ln n})\)的）的时间复杂度，以及 \(O(\sqrt{n})\) 的空间复杂度内求得某些积性函数的前缀和。
一般需要满足 \(f(x)\) 在质数处可以拆成若干个完全积性函数的和，以及可以快速求出 \(f(p^{k})\) 。
一般用来跑 \(n\leq 10^10\) 的数据，这个也可以作为在一些数论题中猜测是否需要用到 Min_25 的依据。
还是从版题引入吧。
简述题意：求积性函数 \(f(x)\) 的前缀和，其中满足 \(f(p^k)=p^k(p^k-1)\) 。
1.分类
首先，把 \(f(x)\) 分为质数和合数，可以得到：
其中 \(LPF\) 是最小质因数（Lowest Prime Factor）的意思。
乍一看似乎没什么意义，为什么要这样子分类需要到后面才能知道。
2.解决质数部分
柿子要挑软的捏
对于上面这个式子，小学生都知道左边的质数部分比较好计算。
可是该怎么算呢？
神仙操作来了！
我们掏出一个dp
定义函数 \(g(n,j)\) ，其满足 \(g(n,j)=\sum^{n}_{i=2}[i\in prim | LFT(i)>p_{j}]f'(x)\) 。
其中的 \(f'(x)\) 应当为一个完全积性函数，并且在质数处的取值与 \(f(x)\) 相同。
注意！，这里的下标从 \(2\) 开始，是不包括 \(i=1\) 的时候的答案的，这是为了方便计算。
为啥要这样定义 \(g\) ？\(f'(x)\) 为什么非得是完全积性函数？一开始我看到这里一头雾水，不过学完之后就豁然开朗了，因此，还请耐心看下去。
注意到 \(f(p^k)=p^{k}(p^{k}-1)=p^{2k}-p^{k}\)，前后两项是可以分开来计算的，而且是一个完全积性函数。由于计算方法一样，下文就只讲一般情况了。
定义完全积性函数 \(f'(x)=x^{k}\) （\(k\) 是任意的常数，在本题中 \(k\) 分别取 \(1\) 和 \(2\) ）。
那么有 \(g(n,j)=\sum^{n}_{i=1}[i\in prim | LFT(i)>p_{j}]x^k\)
接下来试图从 \(g(n,j-1)\) 转移到 \(g(n,j)\) 。
首先，如果有 \(p_{j}>\sqrt{n}\) ，那么显然有 \(g(n,j)=g(n,j-1)\) 。
否则，考虑到 \(g(n,j)\) 应该是比 \(g(n,j-1)\) 少了一些 LFT 为 \(p_{j}\) 的数的贡献，设这个差值为 \(del\) 。
那么可以得到这样的转移式子：
首先可以知道当 \(p_{j}\geq \sqrt{n}\) 的时候是没有意义的，因此第二维的个数不会超过 \(\sqrt{n}\) 。
总状态数是 \(O(n\sqrt{n})\) 级别的，看上去还非常“怪”的样子。
实际上，这里的复杂度上界是 \(O(\frac{n^{\frac{3}{4}}}{\ln n})\)
具体的分析留到后面。
3. 得到总答案
照着上面的思路，在设一个函数 \(S(n,j)\)其和上面的 \(g\) 的定义差不多，不过值由 \(i^k\) 变成了 \(f(i)\)
\(S(n,j)=\sum\limits^{n}_{i=1}[LPF(i)>p_{j}]f(i)\)
在前面已经得到了 \(k=1,2\) 时刻的数组 \(g_{1},g_{2}\)，后面就把这两个作差之后合并为 \(g\) 数组了。
可以发现，当有 \(p_{j}\) 为最大的满足 $p_{j}\leq sqrt{n} $ 的时候， \(g(n,j)\) 中有贡献的值只有质数了。
为了方便，后面的 \(g\) 只写第一项，第二项默认为上面所说的情况。
那么有
答案为 \(S(n,0)\)。
注意到这里有一个 \(e[\not =1]\) 。
因为上面的计算过程中是不统计 \(1\) 的答案的，因此在 \(e \not =1\) 的时候 \(p^{e}\) 这一项会被漏掉，要手动加上去。
而当 \(e=1\) 的时候 \(p^{e}\) 已经在前面的质数部分中被统计过了，因此不需要再加上去。
复杂度分析
引理1:
\(\lfloor \frac{\lfloor\frac{n}{x}\rfloor}{y}\rfloor=\lfloor\frac{n}{xy}\rfloor\)
由于最后我们要求的是 \(S(n,0)\) ，可以发现无论是 \(g\) 数组中还是 \(S\) 数组中第一维出现的数只有可能是 \(\frac{n}{d}\) 的形式，可以知道这样的数只有 \(O(\sqrt{n})\) 个。
但是这些数并不连续，该怎么开数组呢？map
考虑到使用两个数组 \(pos_{1},pos_{2}\)
假设当前需要访问 \(g_{x}\) ，那么对 \(x\) 进行分类。
如果 \(x\leq \sqrt{n}\) ，那么就让它归为 \(pos_{1}\) 中去，对应下标就是 \(x\)。
否则，让它归到 \(pos_{2}\) 中去，对应的下标为 \(\lfloor\frac{n}{x}\rfloor\)。
这一部分的代码实现
ll N;
void init(){
N=sqrt(n);
int &cnt=prim[0];
for(ri i=2;i<=N;++i){
if(!flag[i]) prim[++cnt]=i;
for(ri j=1;j<=cnt&&i*prim[j]<=N;++j){
flag[i*prim[j]]=1;
if(prim[j]%i==0) break;
}
}
for(ll l=1,r;l<=n;l=r+1){//边界情况
r=n/(n/l),w[++tot]=n/r;
g1[tot]=w[tot]%mod;
g2[tot]=dec(g1[tot]*(g1[tot]+1)/2%mod*(2*g1[tot]+1)%mod*inv3%mod,1);//去掉1
g1[tot]=dec(g1[tot]*(g1[tot]+1)/2%mod,1);
if(n/l<=N) pos1[n/r]=tot;//分类
else pos2[n/(n/r)]=tot;
}
}
int getid(ll x){
if(x<=N) return pos1[x];
return pos2[n/x];
}
这样子还有一个好处。
还记得 \(g(n,j)\) 的转移吗？可以发现这样子之后是可以原地滚动的，还可以顺便不用转移第一类情况。
因此保证了时空复杂度的正确性。
引理2:
\([1,n]\) 中的质数个数为 \(O(\frac{n}{\ln n})\) 级别。
第二维中的数为 \([1,\sqrt{n}]\) 中的质数，因此有 \(O(\frac{\sqrt{n}}{\ln \sqrt{n}})\) 个。
为了方便接下来的分析，先给出代码。
Code
#include<cstdio>
#include<cmath>
using namespace std;
#define il inline
#define ri register int
#define ll long long
#define ui unsigned int
il ll read(){
bool f=true;ll x=0;
register char ch=getchar();
while(ch<'0'||ch>'9') {if(ch=='-') f=false;ch=getchar();}
while(ch>='0'&&ch<='9') x=(x<<3)+(x<<1)+(ch^48),ch=getchar();
if(f) return x;
return ~(--x);
}
il void write(const ll &x){if(x>9) write(x/10);putchar(x%10+'0');}
il void print(const ll &x) {x<0?putchar('-'),write(~(x-1)):write(x);putchar('\n');}
il ll max(const ll &a,const ll &b){return a>b?a:b;}
il ll min(const ll &a,const ll &b){return a<b?a:b;}
const ll mod=1e9+7;
const int MAXN=2e5+7;
ll n;
int prim[MAXN];
bool flag[MAXN];
int pos1[MAXN],pos2[MAXN],tot;
ll w[MAXN];
ll g1[MAXN],g2[MAXN];
const ll inv3=333333336;
ll add(ll x,ll y){return (x+=y)<mod?x:x-mod;}
ll dec(ll x,ll y){return (x-=y)<0?x+mod:x;}
ll N;
void init(){
N=sqrt(n);
int &cnt=prim[0];
for(ri i=2;i<=N;++i){
if(!flag[i]) prim[++cnt]=i;
for(ri j=1;j<=cnt&&i*prim[j]<=N;++j){
flag[i*prim[j]]=1;
if(prim[j]%i==0) break;
}
}
for(ll l=1,r;l<=n;l=r+1){//边界情况
r=n/(n/l),w[++tot]=n/r;
g1[tot]=w[tot]%mod;
g2[tot]=dec(g1[tot]*(g1[tot]+1)/2%mod*(2*g1[tot]+1)%mod*inv3%mod,1);//去掉1
g1[tot]=dec(g1[tot]*(g1[tot]+1)/2%mod,1);
if(n/l<=N) pos1[n/r]=tot;//分类
else pos2[n/(n/r)]=tot;
}
}
int getid(ll x){
if(x<=N) return pos1[x];
return pos2[n/x];
}
ll S(ll x,ll y){
if(y&&prim[y]>=x) return 0;
ll r=getid(x),l=getid(y?prim[y]:0),res=dec(dec(g2[r],g1[r]),dec(g2[l],g1[l]));
for(ri i=y+1;i<=prim[0]&&1ll*prim[i]*prim[i]<=x;++i){
for(ll e=1,v=prim[i];v<=x;++e,v=v*prim[i]){
ll now=S(x/v,i);
if(e!=1) now++;
ll t=v%mod;
now=t*(t-1)%mod*now%mod;
res=add(res,now);
}
}
return res;
}
int main(){
n=read();
init();
for(ri i=1;i<=prim[0];++i){
for(ri j=1;j<=tot&&1ll*prim[i]*prim[i]<=w[j];++j){
int nxtl=getid(w[j]/prim[i]),nxtr=(i==1?0:pos1[prim[i-1]]);
g1[j]=dec(g1[j],dec(g1[nxtl],g1[nxtr])*prim[i]%mod);
g2[j]=dec(g2[j],dec(g2[nxtl],g2[nxtr])*prim[i]%mod*prim[i]%mod);
}
}
print(S(n,0)+1);
return 0;
}
该怎么搞我也不清楚
具体复杂度是多少我也不清楚
希腊奶！
数论带师tyy说，处理 \(g\) 的复杂度是 \(O(\frac{n^{\frac{3}{4}}}{\ln n})\) ，而后面求得 \(S\) 部分则是 \(O(n^{1-\epsilon})\)。
反正根据上面的。两个引理也能毛估估出来这玩意是绝对优于线性的
总结
其实到最后看下来，这些操作都是非常行云流水的。
本质上是先利用积性函数的性质先得到了 \(f(x)\) 在质数情况下的答案。
再通过一些巧妙的类似容斥的方法，再根据上面得到的质数答案以及积性函数的性质，快速得到所有的答案。
