杜教筛

杜教筛是一种在非线性时间内求积性函数前缀和的算法。

---

前置知识：

积性函数：

数论函数：定义域为正整数的函数。

积性函数：对于任意互质的整数$a,b$，有$f(ab) = f(a)f(b)$的数论函数

完全积性函数：对于任意整数$a,b$，有$f(ab) = f(a)f(b)$的数论函数

常见的积性函数：$\phi, \mu, \sigma, d$（分别表示：欧拉函数，莫比乌斯函数，因子和函数，因子个数函数）

常见的完全积性函数：$\epsilon, I, id$（$\epsilon(n) = [n = 1], I(n) = 1, id(n) = n$）

狄利克雷卷积：

设$f,g$是两个数论函数，它们的狄利克雷卷积为：$(f*g)(n) = \sum_{d|n}f(d)g(\frac{n}{d})$

运算律：满足交换律，结合律

单位元：$\epsilon$（$f * \epsilon = f$）

几个重要的性质：

1. $\mu * I = \epsilon$
2. $\phi * I = id$
3. $\mu * id = \phi$

利用狄利克雷卷积，我们能很容易的证明莫比乌斯反演：

若：

$$g(n) = \sum_{d|n}f(d) => g = f * I$$

则：

$$f(n) = \sum_{d|n}\mu(d)g(\frac{n}{d}) => \mu * g = \mu * I * f = \epsilon * f = f$$

---

杜教筛：

下面就是正题了，我们现在要求积性函数$f(i)$的前缀和，设$\sum_{i=1}^{n} f(i) = S(n)$

再找一个积性函数$g$，考虑$g * f$的前缀和：

$$\sum_{i=1}^{n} (g * f)(i)$$

$$= \sum_{i=1}^{n} \sum_{d|i} g(d) * f(\frac{i}{d})$$

$$= \sum_{d=1}^{n} g(d) \sum_{i=1}^{\frac{n}{d}} f(\frac{di}{d})$$

$$= \sum_{d=1}^{n} g(d) S(\frac{n}{d})$$

上面这个式子启发我们去考虑$g(1)S(n)$，容易得到：

$$g(1)S(n) = \sum_{i=1}^{n} g(i) S(\frac{n}{i}) - \sum_{j=2}^{n} g(j) S(\frac{n}{j})$$

即：

$$g(1)S(n) = \sum_{i=1}^{n} (g * f)(i) - \sum_{j=2}^{n} g(j) S(\frac{n}{j})$$

至此，我们已经得到了杜教筛的核心式子。

假如我们选取一个合适的积性函数$g$，使得我们能够快速求解$g$的前缀和以及$g * f$的前缀和，那么我们就可以利用数论分块+递归求解上面的式子。

ll GetSum(int n) { // 算 f 前缀和的函数
  ll ans = f_g_sum(n); // 算 f * g 的前缀和
  // 以下这个 for 循环是数论分块
  for(ll l = 2, r; l <= n; l = r + 1) { // 注意从 2 开始
    r = (n / (n / l)); 
    ans -= (g_sum(r) - g_sum(l - 1)) * GetSum(n / l);
    // g_sum 是 g 的前缀和
    // 递归 GetSum 求解
  } return ans; 
}

可以证明，这样做的复杂度是$O(n^{\frac{3}{4}})$的。

我们还可以进一步优化，即先线性筛出前$m$个答案，然后再利用杜教筛求解。

当$m = n^{\frac{2}{3}}$时，复杂度是$O(n^{\frac{2}{3}})$

---

模板题：杜教筛

我们先考虑$\mu$的前缀和，注意到$\mu * I = \epsilon$，于是我们取$g = I, f = \mu$。

那$\epsilon, I$的前缀和都是可以$O(1)$ 求的，直接套用杜教筛的式子即可。

然后是求$\phi$的前缀和，因为有$\phi * I = id$，所以取$g = I, f = \phi$。

这里要注意一下，当$n = 2^{31} - 1$的时候可能会在计算过程中爆$ll$，转换成$ull$运算即可。

#include<bits/stdc++.h>
#include <tr1/unordered_map>
using namespace std;
using namespace std::tr1;

typedef long long ll;
typedef unsigned long long ull;

const int inf = 0x7FFFFFFF;
const int MAXN = 7000000;

int cnt, prim[MAXN + 10];
bool v[MAXN + 10];
ll mu[MAXN], phi[MAXN];
unordered_map<int, ll> ansmu, ansphi;

void init()
{
    v[1] = mu[1] = phi[1] = 1;
    for(int i = 2; i <= MAXN; ++i)
    {
        if(!v[i]) prim[++cnt] = i, mu[i] = -1, phi[i] = i - 1;
        for(int j = 1; j <= cnt && i * prim[j] <= MAXN; ++j)
        {
            v[i * prim[j]] = 1;
            if(i % prim[j]) mu[i * prim[j]] = -mu[i], phi[i * prim[j]] = phi[i] * phi[prim[j]];
            else {mu[i * prim[j]] = 0, phi[i * prim[j]] = phi[i] * prim[j]; break;}
        }
    }
    for(int i = 2; i <= MAXN; ++i) mu[i] += mu[i - 1], phi[i] += phi[i - 1];
}

ll sum_mu(int n)
{
    if(n <= MAXN) return mu[n];
    if(ansmu[n]) return ansmu[n];
    ll ret = 0;
    for(int l = 2, r = 0; r < inf && l <= n; l = r + 1)
    {
        r = n / (n / l);
        ret += (r - l + 1) * sum_mu(n / l);
    }
    return ansmu[n] = 1ll - ret;
}

ll sum_phi(int n)
{
    if(n <= MAXN) return phi[n];
    if(ansphi[n]) return ansphi[n];
    ll ret = 0;
    for(int l = 2, r = 0; r < inf && l <= n; l = r + 1)
    {
        r = n / (n / l);
        ret += (r - l + 1) * sum_phi(n / l);
    }
    return ansphi[n] = (ull)n * (n + 1ull) / 2ull - ret;
}

int main()
{
    init(); int T, n;
    scanf("%d", &T);
    while(T--)
    {
        scanf("%d", &n);
        printf("%lld %lld\n", sum_phi(n), sum_mu(n));
    }
    return 0;
}