设 $f$ 为积性函数,且其在质数处的取值 $f(p)$ 是关于 $p$ 的多项式,且其在质数的幂的取值 $f(p^c)$ 能快速计算,则 Min_25 筛可以在 $\mathcal{O}(\frac{n^{\frac{3}{4}}}{\log n})$ 的时间复杂度内算出 $\sum_{i=1}^nf(i)$。


思路

为了方便,设 $\mathrm{P}$ 为质数集,$p_i$ 为第 $i$ 小的质数。

第一步

先考虑如何求 $\sum_{i=1}^n[i\in\mathrm{P}]f(i)$。

因为 $f(p)$ 是关于 $p$ 的多项式,因此我们考虑对于若干个函数 $f_k(i)=i^k$,计算 $\sum_{i=1}^x[i\in\mathrm{P}]f_k(i)$,然后再加减一下就可以得到 $\sum_{i=1}^x[i\in\mathrm{P}]f(i)$ 了。

考虑怎么求 $\sum_{i=1}^x[i\in\mathrm{P}]f_k(i)$。为了方便,本节剩余部分中直接将 $f_k(i)$ 写作 $f(i)$。

考虑一个 DP。设

$$ g(n,j)=\sum_{i=1}^n[i\in\mathrm{P}\lor\min_{p|i,p\in\mathrm{P}}\left\{p\right\}>p_j]f(i) $$

可以发现,$\sum_{i=1}^n[i\in\mathrm{P}]f(i)$ 就等于 $g(n,|\mathrm{P}|)$。

考虑 $g(n,j)$ 的实际意义。可以发现 $g(n,j)$ 就等于埃筛时用 $j$ 个质数筛后未被筛去的所有数的 $f$ 值之和加上所有质数的 $f$ 值之和。

根据这个实际意义考虑转移。

  • 当 $p_j^2>n$ 时,显然不会筛去任何数,所以此时 $g(n,j)=g(n,j-1)$。
  • 当 $p_j^2\leq n$ 时,筛掉了最小质因子为 $p_j$ 的那些数。考虑将这些数除掉一个 $p_j$,那么除完之后的数最小质因子要 $\geq p_j$,则这些数的 $f$ 值之和为 $g(\frac{n}{p_j},j-1)-\sum_{i=1}^{j-1}f(p_i)$。又因为 $f$ 是完全积性函数,所以此时 $g(n,j)=g(n,j-1)-f(p_j)\times\left[g(\frac{n}{p_j},j-1)-\sum_{i=1}^{j-1}f(p_i)\right]$。

综上可以得到总转移

$$ g(n,j)=\begin{cases}g(n,j-1),&p_j^2>n\\g(n,j-1)-f(p_j)\times\left[g(\frac{n}{p_j},j-1)-\sum_{i=1}^{j-1}f(p_i)\right],&p_j^2\leq n\end{cases} $$

边界为 $g(n,0)=\sum_{i=1}^nf(i)$。

第二步

$$ S(n,j)=\sum_{i=1}^n\left[\min_{p|i,p\in\mathrm{P}}\left\{p\right\}\geq p_j\right]f(i) $$

考虑转移,分质数和合数两部分计算,可以得到

$$ S(n,j)=g(n,|\mathrm{P}|)-\sum_{i=1}^{j-1}f(p_i)+\sum_{k\geq j,p_k^2\leq n}\sum_{e\geq 1,p_k^{e+1}\leq n}\left[S\left(\frac{n}{p_k^e},k+1\right)\times f(p_k^e)+f(p_k^{e+1})\right] $$

合数部分的意思是,枚举最小质因数以及次数,由于是积性函数所以可以直接乘,然后少算了 $f(p_k^2),f(p_k^3),\cdots$ 所以需要补进来($f(p_k)$ 在质数部分计算过了)。

那么最后答案即为 $S(n,1)+f(1)$。

细节

递归求解 $S(n,j)$。则递归到的那些 $n$ 为 $\left\lfloor\frac{n}{a}\right\rfloor,\left\lfloor\frac{\left\lfloor\frac{n}{a}\right\rfloor}{b}\right\rfloor,\cdots$。然后有一个这样的定理

$$ \left\lfloor\frac{\left\lfloor\frac{n}{a}\right\rfloor}{b}\right\rfloor=\left\lfloor\frac{n}{ab}\right\rfloor $$

因此递归到的 $n$ 为 $\left\lfloor\frac{n}{i}\right\rfloor$,而这样的值的数量是 $\mathcal{O}(\sqrt{n})$ 级别的。

于是我们仅需要计算出 $\mathcal{O}(\sqrt{n})$ 个 $g(n,|\mathrm{P}|)$ 即可。

计算 $g(n,|\mathrm{P}|)$ 时,当 $p_j>n$ 时 $g$ 值不会再改变,因此只需要预处理出 $\sqrt{n}$ 内的素数即可。

实现

首先筛出 $\sqrt{n}$ 内的质数。

然后求出所有 $\left\lfloor\frac{n}{i}\right\rfloor$ 的值并存下来。存的时候开两个数组,如果 $\left\lfloor\frac{n}{i}\right\rfloor\leq\sqrt{n}$ 则以 $\left\lfloor\frac{n}{i}\right\rfloor$ 为下标存到第一个数组中,否则以 $\left\lfloor\frac{n}{\left\lfloor\frac{n}{i}\right\rfloor}\right\rfloor$ 为下标存到第二个数组中。

接着对求出的这些值,求出它对应的 $g(n,|\mathrm{P}|)$。

最后直接递归求解 $S(n,1)$ 即可,记得加上 $f(1)$。

理解

本部分不一定准确,如有错误烦请指正。

由于最后算出来的 $g(n,|\mathrm{P}|)$ 仅包含质数处的取值,所以我们假装 $f$ 在非质数处的取值和质数处相同是没有影响的。

在计算 $g(n,|\mathrm{P}|)$ 时的转移要求 $f$ 是完全积性函数,所以需要拆成若干个函数分开求。这就是要求 $f(p)$ 为关于 $p$ 的多项式的原因。

计算 $S(n,j)$ 时的转移要求 $f$ 是积性函数且 $f(p^c)$ 能够快速计算。

例子

【模板】Min_25 筛

$f(p)=p^2-p$,因此分别对 $f_1(i)=i$ 和 $f_2(i)=i^2$ 求 $g(n,|\mathrm{P}|)$,答案为 $S(n,1)+1$。

// ===================================
//   author: M_sea
//   website: http://m-sea-blog.com/
// ===================================
#include <bits/stdc++.h>
#define re register
using namespace std;
typedef long long ll;

inline ll read() {
    ll X=0,w=1; char c=getchar();
    while (c<'0'||c>'9') { if (c=='-') w=-1; c=getchar(); }
    while (c>='0'&&c<='9') X=X*10+c-'0',c=getchar();
    return X*w;
}

const int N=200000+10;
const int mod=1e9+7,inv6=166666668,inv2=500000004;

ll n; int s;

int p[N],v[N],cnt=0;
inline void sieve(int n) {
    for (re int i=2;i<=n;++i) {
        if (!v[i]) p[++cnt]=i;
        for (re int j=1;j<=cnt&&i*p[j]<=n;++j) {
            v[i*p[j]]=1;
            if (i%p[j]==0) break;
        }
    }
}

ll w[N]; int id1[N],id2[N]; int top=0;
int gs1[N],gs2[N],g1[N],g2[N];

inline int S(ll x,int j) {
    if (x<=1||p[j]>x) return 0;
    int k=x<=s?id1[x]:id2[n/x];
    int res=((1ll*g2[k]-g1[k]-gs2[j-1]+gs1[j-1])%mod+mod)%mod;
    for (re int i=j;i<=cnt&&1ll*p[i]*p[i]<=x;++i) {
        ll p1=p[i],p2=1ll*p[i]*p[i];
        for (re int e=1;p2<=x;++e,p1=p2,p2*=p[i]) {
            int q1=p1%mod,q2=p2%mod;
            res=(res+1ll*S(x/p1,i+1)*q1%mod*(q1-1)+1ll*q2*(q2-1))%mod;
        }
    }
    return res;
}

int main() {
    n=read(),s=sqrt(n); sieve(s);
    for (re int i=1;i<=cnt;++i) {
        gs1[i]=(gs1[i-1]+p[i])%mod;
        gs2[i]=(gs2[i-1]+1ll*p[i]*p[i])%mod;
    }
    for (re ll l=1,r;l<=n;l=r+1) {
        r=n/(n/l); w[++top]=n/l;
        n/l<=s?id1[n/l]=top:id2[n/(n/l)]=top;
        int x=(n/l)%mod;
        g1[top]=(1ll*x*(x+1)%mod*inv2%mod-1+mod)%mod;
        g2[top]=(1ll*x*(x+1)%mod*(x+x+1)%mod*inv6%mod-1+mod)%mod;
    }
    for (re int i=1;i<=cnt;++i)
        for (re int j=1;j<=top&&1ll*p[i]*p[i]<=w[j];++j) {
            ll x=w[j]/p[i]; int k=x<=s?id1[x]:id2[n/x];
            g1[j]=(g1[j]-1ll*p[i]*(g1[k]-gs1[i-1]+mod)%mod+mod)%mod;
            g2[j]=(g2[j]-1ll*p[i]*p[i]%mod*(g2[k]-gs2[i-1]+mod)%mod+mod)%mod;
        }
    printf("%d\n",(S(n,1)+1)%mod);
    return 0;
}

版权属于:M_sea

本文链接:https://m-sea-blog.com/archives/min_25/

转载时须注明出处及本声明

最后修改:2020 年 02 月 14 日 11 : 06 AM