扯淡

min_25篩是由min_25提出的求積性函數前綴和的亞線性算法，和一個叫“擴展埃氏篩”的東西有著微妙的關系。

至于是什么關系，我也不太清楚，反正有人說很像有人說就是一個東西（霧）

這段話并不是廢話

約定

為了方便后面描述，這里寫一些用到的約定和符號表示，以免產生恐懼

$1$被開除正整數籍 也就是說“前綴和”之類的都是從$2$開始，對答案所求的前綴和同樣，最后手動加$1$

$\pi (n)$表示$1～n$中質數個數的規模（其實問題不大，后面就懂了）

$p_i$表示第$i$個質數，單獨的$p$均表示質數

$prime$表示質數集合 $minp(i)$表示$i$的最小質因子

$p^c$表示一個只含一個質因子的數

流程

首先所求的函數$f(n)$需要滿足：

是個積性函數

在質數處的取值$f(p)$是一個關于$p$的多項式。我們把這個多項式拆成若干單項式分別計算在相加，這樣就變成了"$f(p)$的值是一個關于$p$的單項式"，我們記為$f(p)=p^k$

$f(p^c)$的值可以快速計算

首先考慮計算這個東西

$$\sum _{i=2}^n[i\in prime]i^k$$

注意是$i^k$不是$f(i)$，雖然這里還沒有區別，但它們只是質數位置相等

首先既然叫擴展埃氏篩先回憶一下埃氏篩在干什么：開始寫下所有正整數，然后從小到大，如果一個沒被篩說明它是質數，用它篩掉后面的倍數

我們可以用這個思路計算上面的式子，先算出所有的和，再用質數篩掉所有合數項

設$g(n,j)$表示用前$j$個質數篩了之后剩余項的和

$$g(n,j)=\sum _{i=2}^n[i\in primeorminp(i)>p_j]i^k$$

和就是$g(n,\pi (\sqrt{n}))$

注意合數的最小質因子不會超過$\sqrt{n}$,所以如果$p_j^2>n$有$g(n,j)=g(n,j?1)$

對剩下情況考慮轉移，即從$g(n,j?1)$去掉被$p_j$刪掉的項

$$g(n,j?1)=\sum _{i=2}^n[i\in primeorminp(i)\ge p_j]$$

如果腦子轉不過來，可以強行算

$[i\in primeorminp(i)\ge p_j]$且不滿足$[i\in primeorminp(i)>p_j]$
即
$[i\in primeorminp(i)\ge p_j]$且$[i\in /prime]且[minp(i)\le p_j]$
綜上
$[i\in /prime,minp(i)=p_j]$

所以

$$g(n,j)=g(n,j?1)?\sum _{i=2}^n[i\in /prime,minp(i)=p_j]i^k$$

可以提一個$p_j$出來,因為不能有質數，就把$p_j$去掉，剛好是$1$

$$g(n,j)=g(n,j?1)?p_j^k\sum _{i=2}^{?\frac{n}{p_j}?}[minp(i)\ge p_j]i^k$$

再次觀察

$$g(n,j)=\sum _{i=2}^n[i\in primeorminp(i)>p_j]i^k$$

發現可以拆成小于$p_j$的質數和大于等于$p_j$的所有數，第二個就是剛才的式子

$$g(n,j)=g(n,j?1)?p_j^k(g(?\frac{n}{p_j}?,j?1)?\sum _{i=1}^{j?1}{p}_i^k)$$

注意$n$由于都是一直整除，所以只會有$O(\sqrt{n})$種取值，可以整除分塊找出來強行離散化。在記錄某個值$v$所在位置的時候，如果$v>\sqrt{n}$,我們另開一個數組存到$?\frac{n}{v}?$里面

然后后面只會用到最后一項，所以第二維可以滾掉

求答案時，設

$$S(n,j)=\sum _{i=2}^n[minp(i)>p_j]f(i)$$

分質數和合數分別計算

質數部分用$g$算出來的減去前$j$個

$$g(n,\pi (\sqrt{n}))?\sum _{i=1}^j{p}_i^k$$

合數部分枚舉最小質因子和它的次數

$$\sum _{k=j+1}^{p_k\le \sqrt{n}}\sum _{e=1}^{p_k^e\le n}f(p_k^e)(S(?\frac{n}{p_k^e}?,k)+[e>1])$$

$[e>1]$指如果指數大于$1$它本身就是合數，需要統計答案

兩部分相加即可

總復雜度大概$O(n^{\frac{2}{3}})$常數較小，大約$2s$過$1e10$,$3s$過$1e11$

模板題

#include <iostream> #include <cstdio> #include <cstring> #include <cctype> #include <cmath>#define MAXN 200005 using namespace std; typedef long long ll; const int MOD=1e9+7,INV6=(MOD+1)/6; inline int add(const int& x,const int& y){return x+y>=MOD? x+y-MOD:x+y;} inline int dec(const int& x,const int& y){return x<y? x-y+MOD:x-y;} int np[MAXN],pl[MAXN],cnt; inline void init(const int& N) {np[1]=1;for (int i=2;i<=N;i++){if (!np[i]) pl[++cnt]=i;int x;for (int j=1;(x=i*pl[j])<=N;j++){np[x]=1;if (i%pl[j]==0) break;}} } ll val[MAXN]; int tot; int g1[MAXN],g2[MAXN],sum1[MAXN],sum2[MAXN]; int key[MAXN],yek[MAXN]; ll n; int m; inline int getkey(const ll& v){return v<=m? key[v]:yek[n/v];} int S(ll n,int j) {if (pl[j]>=n) return 0;int k=getkey(n);int ans=dec(dec(g2[k],g1[k]),dec(sum2[j],sum1[j]));for (int k=j+1;k<=cnt&&(ll)pl[k]*pl[k]<=n;k++)for (ll e=1,pe=pl[k];pe<=n;e++,pe*=pl[k])ans=add(ans,pe%MOD*(pe%MOD-1)%MOD*(S(n/pe,k)+(e>1))%MOD); return ans; } int main() {scanf("%lld",&n);m=sqrt(n);init(m);for (ll l=1,r;l<=n;l=r+1){r=n/(n/l);val[++tot]=n/l;int t=val[tot]%MOD;g1[tot]=((ll)t*(t+1)/2)%MOD-1;g2[tot]=(ll)t*(t+1)%MOD*(2*t+1)%MOD*INV6%MOD-1;if (val[tot]<=m) key[val[tot]]=tot;else yek[n/(n/l)]=tot;}for (int j=1;j<=cnt;j++){for (int i=1;i<=tot&&(ll)pl[j]*pl[j]<=val[i];i++){int k=getkey(val[i]/pl[j]);g1[i]=dec(g1[i],(ll)pl[j]*dec(g1[k],sum1[j-1])%MOD);g2[i]=dec(g2[i],(ll)pl[j]*pl[j]%MOD*dec(g2[k],sum2[j-1])%MOD);} sum1[j]=add(sum1[j-1],pl[j]),sum2[j]=add(sum2[j-1],(ll)pl[j]*pl[j]%MOD);}printf("%d\n",S(n,0)+1);return 0; }

總結

以上是生活随笔為你收集整理的min_25筛详解的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。