杜教筛

前置知识

杜教筛

求积性函数 \(f(n)\) 的前缀和，即 \(\sum\limits_{i = 1}^{n} f(i)\)。

考虑狄利克雷卷积能否帮助我们优化。

\[\begin{aligned} \sum_{i = 1}^{n}(f * g)(i) &= \sum_{i = 1}^{n}\sum_{d \mid i}f(d)g\left(\frac{i}{d}\right) \\ &= \sum_{d = 1}^{n}\sum_{i = 1}^{\lfloor\frac{n}{d}\rfloor}g(d)f(i) \\ &= \sum_{d = 1}^{n}g(d)\left(\sum_{i = 1}^{\lfloor\frac{n}{d}\rfloor}f(i)\right) \end{aligned} \]

记 \(S(n) = \sum\limits_{i = 1}^{n}f(i)\)，则有：

\[\sum_{i = 1}^{n}(f * g)(i) = \sum_{d = 1}^{n}g(d)S\left(\left\lfloor\frac{n}{d}\right\rfloor\right) \]

欲求 \(S(n)\)，可令 \(d = 1\) 发现：

\[g(1)S(n) = \sum_{i = 1}^{n}(f * g)(i) - \sum_{d = 2}^{n}g(d)S\left(\left\lfloor\frac{n}{d}\right\rfloor\right) \]

如果我们能找到一个函数 \(g(n)\) 使得 \((f * g)(n)\) 与 \(g(n)\) 的前缀和能够快速求出，那么 \(f(n)\) 的前缀和就可以通过数论分块快速求出。

ll get_sum(ll n)
{
    if (n < N) return pref[n]; // 预处理的前缀和
    if (sumf[n]) return sumf[n]; // 哈希存储前缀和
    ll res = get_sum(n); // f * g 的前缀和
    for (ll l = 2, r; l <= n; l = r + 1)
    {
        r = n / (n / l);
        res -= get_g(l, r) * get_sum(n / l); // g(l ~ r) * f(n / l)
    }
    return sumf[n] = res;
}

\(\varphi(n)\) 的前缀和

根据 \(\varphi * {\mathbf{1}} = {\rm{Id}}\)，容易得出以下结论：

\[\sum_{i = 1}^{n} {\mathbf{1}} = n \]

\[\sum_{i = 1}^{n} {\rm{Id}}(n) = \frac{n(n+1)}{2} \]

用哈希记忆化一下，考虑递归求解，递归代码如下：

ll get_sumphi(ll n)
{
    if (n < N) return phi[n];
    if (sumphi[n]) return sumphi[n];
    ll res = n * (n + 1) / 2; // f * g 的前缀和
    for (ll l = 2, r; l <= n; l = r + 1)
    {
        r = n / (n / l);
        res -= (r - l + 1) * get_sumphi(n / l); // g(l ~ r) * f(n / l)
    }
    return sumphi[n] = res;
}

\(\mu(n)\) 的前缀和

根据 \(\mu * {\mathbf{1}} = \varepsilon\)，容易得出以下结论：

\[\sum_{i = 1}^{n} {\mathbf{1}} = n \]

\[\sum_{i = 1}^{n} \varepsilon(n) = 1 \]

同样的，可以写出下面的代码：

ll get_summu(ll n)
{
    if (n < N) return mu[n];
    if (summu[n]) return summu[n];
    ll res = 1;
    for (ll l = 2, r; l <= n; l = r + 1)
    {
        r = n / (n / l);
        res -= (r - l + 1) * get_summu(n / l);
    }
    return summu[n] = res;
}

题目实例

洛谷 P4213 【模版】杜教筛

给定一个正整数，求

\[ans_1=\sum_{i=1}^n\varphi(i) \]

\[ans_2=\sum_{i=1}^n \mu(i) \]

对于全部的测试点，保证 \(1 \leq T \leq 10\)，\(1 \leq n \lt 2^{31}\)。

也就是我们上面讨论的 \(\varphi(n)\) 与 \(\mu(n)\) 的前缀和。

参考代码

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = 1e7 + 10;
int cnt, st[N], primes[N];
ll n, phi[N], mu[N];
unordered_map<ll, ll> sumphi, summu;

ll get_sumphi(ll n)
{
    if (n < N) return phi[n];
    if (sumphi[n]) return sumphi[n];
    ll res = n * (n + 1) / 2;
    for (ll l = 2, r; l <= n; l = r + 1)
    {
        r = n / (n / l);
        res -= (r - l + 1) * get_sumphi(n / l);
    }
    return sumphi[n] = res;
}

ll get_summu(ll n)
{
    if (n < N) return mu[n];
    if (summu[n]) return summu[n];
    ll res = 1;
    for (ll l = 2, r; l <= n; l = r + 1)
    {
        r = n / (n / l);
        res -= (r - l + 1) * get_summu(n / l);
    }
    return summu[n] = res;
}

void init()
{
    phi[1] = mu[1] = 1;
    for (ll i = 2; i < N; i ++ )
    {
        if (!st[i]) primes[cnt ++ ] = i, phi[i] = i - 1, mu[i] = -1;
        for (ll j = 0; primes[j] * i < N && j < cnt; j ++ )
        {
            st[primes[j] * i] = 1;
            if (i % primes[j] == 0)
            {
                phi[primes[j] * i] = phi[i] * primes[j];
                break;
            }
            phi[primes[j] * i] = phi[i] * (primes[j] - 1);
            mu[primes[j] * i] = -mu[i];
        }
    }
    for (int i = 1; i < N; i ++ ) phi[i] += phi[i - 1], mu[i] += mu[i - 1];
}

void solve()
{
    scanf("%lld", &n);
    printf("%lld %lld\n", get_sumphi(n), get_summu(n));
}

int main()
{
    int T;
    init();
    scanf("%d", &T);
    while (T -- ) solve();
    return 0;
}

洛谷 P3768 简单的数学题

\[\left(\sum_{i=1}^n\sum_{j=1}^n ij \gcd(i,j)\right) \bmod p \]

对于所有的数据，\(n \leq 10^{10}\)。

对于所有的数据，\(5 \times 10^8 \leq p \leq 1.1 \times 10^9\) 且 \(p\) 为质数。

对 \(\gcd{(i, j)}\) 进行欧拉反演：

\[\begin{aligned} & \sum_{i=1}^n\sum_{j=1}^n ij \gcd(i,j) \\ =& \sum_{i=1}^n\sum_{j=1}^n ij \sum_{d \mid \gcd{(i, j)}} \varphi(d) \\ =& \sum_{i=1}^n\sum_{j=1}^n ij \sum_{d \mid i, d \mid j} \varphi(d) \\ =& \sum_{d=1}^n\varphi(d) \sum_{d \mid i}^n\sum_{d \mid j}^n ij \\ =& \sum_{d=1}^n\varphi(d)d^2 \sum_{i=1}^{\lfloor\frac{n}{d}\rfloor}\sum_{j=1}^{\lfloor\frac{n}{d}\rfloor} ij \\ =& \sum_{d=1}^n\varphi(d)d^2 \left(\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor} i\right)^2 \\ =& \sum_{d=1}^n\varphi(d)d^2 \left(\sum_{i=1}^{\lfloor\frac{n}{d}\rfloor} i\right)^2 \\ \end{aligned} \]

记 \(f(n) = n^2\varphi(n)\)，\(g(n) = n^2\)。

\[\begin{aligned} (f * g)(n) &= \sum_{d \mid n}d^2\varphi(d)\left(\frac{n}{d}\right)^2 \\ &= n^2\sum_{d \mid n}\varphi(d) \\ &= n^3 \end{aligned} \]

记 \(S(n) = \sum\limits_{i=1}^{n}i^2\varphi(i)\)，\(h(n) = \left(\sum\limits_{i=1}^n i\right)^2 = \frac{n^2(n+1)^2}{4}\)，对后面这个式子进行数论分块，杜教筛求前面式子即可，同样的预处理出 \(S(n)\) 的前 \(n^{\frac{2}{3}}\) 项可以做到 \(O(n^{\frac{2}{3}})\)。