min25筛

Jul 30, 2019 由 小羊

这也是个用于求积性函数前缀和的筛法。

对于积性函数 $f(n)$，满足这样几个要求时可以用min25来做。

• 当 $n \in prime$ 时， $f(n)$ 可以使用关于 $n$ 的多项式表示。

• 当 $n \in prime$ 时，$f(n^k)$ 要能够快速计算答案。

预处理 Σf(素数)

在min25筛的过程中，我们需要预处理多项式的 $x^k$ 相关的这个东西

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

那么记所有质数的集合为 $P$，第 $j$ 个素数为 $P_j$，$x$ 的最小质因子为 $minp(x)$，我们需要预处理这样一个东西

$g(n,j) = \sum_{i=2}^n i^k [i \in P ~ or ~ minp(i) > P_j]$

举个例子吧

$\begin{aligned} g(n,0) &= 2^k + 3^k + 4^k + 5^k + 6^k + 7^k + \cdots + n^k \\ g(2n,1) &= 2^k + 3^k + 5^k + 7^k + 9^k + 11^k + \cdots + (2n-1)^k \\ g(6n,2) &= 2^k + 3^k + 5^k + 7^k + 11^k + \cdots + (6n-5)^k + (6n-1)^k \end{aligned}$

相当于在这个求和中，$g(n,j)$ 将前 $j$ 个素数的倍数全部去掉了，只留下了那些素数的 $k$ 次方和。

$\begin{aligned} g(2n,0) - g(2n,1) &= 4^k + 6^k + 8^k + 10^k + 12^k + \cdots + (2n)^k \\ &= 2^k \left( 2^k + 3^k + 4^k + 5^k + 6^k + \cdots + n^k \right) \\ \\ g(6n,1) - g(6n,2) &= 9^k + 15^k + 21^k + 27^k + 33^k + \cdots + (6n-3)^k \\ &= 3^k \left( 3^k + 5^k + 7^k + 9^k + 11^k + \cdots + (2n-1)^k \right) \end{aligned}$

那么从 $g(n,j-1)$ 转移到 $g(n,j)$ 的过程中，我们考虑去掉 $p_j^k$ 的倍数 $\sum (p_j \prod p_{j+i}^{e_i})^k$。

差不多就是有这样一个转移

$g(n,j) = g(n,j-1) - p_j^k \left( g\left(\left\lfloor\frac{n}{p_j}\right\rfloor, j-1\right) - \sum_{x=1}^{j-1} p_x^k \right)$

为什么要减去前 $j-1$ 个素数 $k$ 次方和的转移呢？因为 $(p_x p_j)^k$ 在这之前都被筛掉了啊。

如果 $p_j^2 > n$，那肯定没有东西能被删掉了。所以有这样一个转移

$g(n,j) = \begin{cases} g(n,j-1), &p_j^2 \gt n \\ g(n,j-1) - p_j^k \left( g(\lfloor n / p_j\rfloor, j-1) - \sum_{x=1}^{j-1} p_x^k \right), &p_j^2 \le n \end{cases}$

我们需要一些转换，才能把上面那个东西从递归变成迭代。

我们会发现，我们可能需要经常性的计算 $g(\lfloor n/x\rfloor, j)$ 的值，那我们就一起维护了吧。

这里是预处理的代码。

// euler sieve, spk for sum p^k from p_1 ~ p_pcn
int pri[MAXN], pcn, spk[MAXN];
int m, Sqrt, id1[MAXN], id2[MAXN];
long long gk[MAXN], w[MAXN];

inline long long sum_ik(lld n); // for i=2..n, sum i^k

void init()
{
    m = 0, Sqrt = sqrt(n);

    for (long long L = 1, R; L <= n; L = R+1)
    {
        R = n / (n / L), w[++m] = n / L;
        (w[m] <= Sqrt ? id1[w[m]] : id2[R]) = m;
        gk[m] = sum_ik(w[m]);
    }

    for (int j = 1; j <= pcn; j++) // p_pcn^2 should be greater than n
    {
        for (int i = 1; i <= m && pri[j] <= w[i]/pri[j]; i++)
        {
            int d = (w[i]/pri[j]<=Sqrt) ? id1[w[i]/pri[j]] : id2[n/(w[i]/pri[j])];
            gk[i] -= (spk[j] - spk[j-1]) * (gk[d] - spk[j-1]);
        }
    }
}

其中，id1[i] 表示数论分块中整除值 $\lfloor n/l \rfloor = i$ 的块编号，id2[n/i] 表示整除值 $\lfloor n / l \rfloor = i$ 的块编号，w[s] 表示数论分块第 $s$ 个块的整除值。这样比较节省空间。

然后考虑下面那个迭代。首先，$w_i$ 就是 g 函数的第一个参数。

当 $p_j^2 > w_i$，也就是第 $i$ 个块整除值的时候，没有必要继续处理了。内层的 $i=1 \cdots m$ 相当于从大往小枚举块，也就不用利用滚动数组什么的记录上一轮 $j-1$ 时的具体值了。

最后，init执行完时，gk[i] 就是第 $i$ 个整除值 $w_i$ 时，$g(w_i, pcn)$ 的值了，也就是我们想要的东西了。

由于那个多项式可能会有好几个 $p^x$，所以都维护出来吧。

求 Σf(i)

我们刚才去掉了合数处的点值，那么我们还要根据积性函数的性质把那些点值加回来。

那么我们现在构造这样一个函数

$S(n,j) = \sum_{x=2}^n f(x)~[ minp(x) \ge p_j ]$

例如

$\begin{aligned} S(20,3) &= f(5) + f(7) + f(11) + f(13) + f(17) + f(19) \\ S(20,2) - S(20,3) &= f(3) + f(9) + f(15) \\ S(20,1) - S(20,2) &= f(2) + f(4) + f(6) + \cdots + f(18) + f(20) \end{aligned}$

我们需要从最大的质数开始，一个个往回补上点值。

$S(n,j) = \sum_{k \ge j} \sum_{e} f(p_k^e)\left(1 + S\left(\left\lfloor n/p_k^e\right\rfloor, k+1\right)\right)$

考虑 $p_k^e$ 的倍数，显然当 $e=0$ 的时候不能算上贡献；当 $e>0$ 的时候，$\sum p_k^e \prod p_{k+i}^{e_i}$ 的贡献由 $f(p_k^e)$ 乘上 $\sum f(\prod p_{k+i}^{e_i})$ 得到。为什么还有个 $1+$ 呢？$S(x,k)$ 里面可没有 $p_j$ 的东西，不加就算不上 $\sum f(p_k^e)$ 的贡献了呀。

那我们整理一下这个鬼一样的式子。

$\begin{aligned} S(n, j) &= \sum_{k \ge j} f(p_k) + \sum_{k \ge j} \sum_{e} \left( f(p_k^e) S(\lfloor n/p_k^e\rfloor, k+1) + f(p_k^{e+1}) \right) \\ &= g(n, |P|) - \sum_{i=1}^{j-1} f(p_i) + \sum_{k \ge j} \sum_{e} \left( f(p_k^e) S(\lfloor n/p_k^e\rfloor) + f(p_k^{e+1}) \right) \end{aligned}$

所以这一部分的代码可以写成

long long S(long long x, int y)
{
    if (x <= 1 || pri[y] > x) return 0;
    int d = x <= Sqrt ? id1[x] : id2[n/x];
    long long ans = g(x) - g(pri[y]-1);

    for (int i = y; i <= pcn && pri[i] <= x / pri[i]; i++)
    {
        long long p1 = pri[i], p2 = p1 * p1;
        for (int e = 1; p1 <= x / pri[i]; ++e, p1 = p2, p2 *= pri[i])
            ans += S(x/p1, i+1) * f(i,e) + f(i,e+1);
    }

    return ans;
}

其中 g(n) 是 $\sum_{i=1}^n f(i) [i \in prime]$，也就是上面式子中的 $g(n,|P|)$ 可以从上面预处理的式子中得到。而 f(i,e) 就是 $f(p_i^e)$。

最后你要求的 $\sum_{i=1}^n f(i) = S(n,1) + 1$

举例来说，我已知积性函数 $f(p^e) = C_{e+A-1}^{A-1}$，那么 $f(p)=A$，我可以预处理上面预处理部分 sp0, g0 等内容，那么 g(x) - g(pri[y]-1) 应该写为 A * (g0[d]-(y-1))，而 f(i,e) 则写为 C(e+A-1, A-1)。对应例题为 HDU6537。

举例来说，我已知积性函数 $f(p^e) = Ae+1$，那么 $f(p)=A+1$，我可以预处理 sp0, g0 等内容，那么 g(x) - g(pri[y]-1) 应该写为 (A+1) * (g0[d]-(y-1))，而 f(i,e) 则写为 A*e+1。对应具体例题为 SPOJ DIVCNT3。

再举例来说，我已知积性函数 $f(p^e) = p \oplus e$，那么可以在 $p=2$ 时 $f(2)=3$，然后其他时候有 $f(p)=p-1$，预处理 sp0, sp1, g0, g1 等内容，那么 g(x) - g(pri[y]-1) 应该写为 g1[d]-g0[d]-sp1[y-1]+y-1，并特判 if (y == 1) ans += 2;，而 f(i,e) 则写为 pri[i]^e。对应具体例题为 LOJ6053。

后记

所以说这是不是很像埃式筛啊，所以它又被称为拓展埃式筛。

据说时间复杂度是

$O\left( \frac{n^{\frac{3}{4}}}{\log n} \right)$

咱也不会证，咱也不敢证。