杜教筛 学习笔记

杜教筛用于计算某一类积性函数的前缀和。

思想

构造 Dirchlet 卷积，把一个困难的前缀和转化为两个相对简单的前缀和。

形式化地说，我们有数论函数 $f,g$，记 $F(n) = \sum\limits_{i = 1}^n f(i)$ ，那么有：
$$
\begin{aligned}
\sum\limits_{i = 1}^n fg(i) &= \sum\limits_{i=1}^n \sum_{x\times y = i}f(x)\times g(y) \
&= \sum\limits_{y=1}^n g(y)\sum\limits_{x\times y \leq n} f(x)\
&= \sum\limits_{y=1}^n g(y)\times F(\lfloor \frac{n}{y}\rfloor)\
\end{aligned}
$$
也就是说，我们有：
$$
g(1)F(n) = \sum\limits_{i = 1}^n fg(i) - \sum\limits_{y=2}^n g(y)\times F(\lfloor \frac{n}{y}\rfloor)
$$
这里如果我们假设 $f g$ 的前缀和和 $g$ 的区间和都可以 $O(1)$ 计算（这也是我们应用杜教筛的一个条件），那么计算 $F$，我们需要的就是计算后面 $\sum\limits_{y=2}^n g(y)\times F(\lfloor \frac{n}{y}\rfloor)$ ，那么这个显然可以数论分块，如果直接数论分块，我们*递归计算 $F(\lfloor \frac{n}{y}\rfloor) $就是 $\mathcal{O}(\sqrt n)$ 次的了。

好，我们解决了单次计算的循环次数问题，那么我们发现每次递归的 $F(x)$ 的取值实质上是 $1,2\dots \sqrt n$ ，这个如果直接记忆化暴力计算，通过调和级数我们可以知道复杂度是 $\mathcal{O}(n^\frac{3}{4})$ 的。

但是我们发现，对于较小的 $n$ 通过 $\mathcal{O}\sqrt n$ 来计算是血亏的，我们运用一个类似于根号分治的思想，对于 $1\dots n^c$ 的答案用 $\mathcal{O}(n^c)$ 的时间预处理出来，通过一些积分魔术我们可以得到当 $c = \frac{2}{3}$ 的时候复杂度最优秀，为 $\mathcal{O}(n^{ \frac{2}{3}})$ 。

记忆化的部分，我们可以直接记录一个 ans[x] 代表 $\frac{n}{x}$ 时候的答案。对于计算 $f *g$ 以及 $g$ 的前缀和，我们发现只需要在 $\mathcal{O}(\sqrt n)$ 的时间内算完就不影响复杂度，实际上你也可以套一层杜教筛（

例题

P4213 【模板】杜教筛（Sum)

给定一个正整数，求

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

$$ans_2=\sum_{i=1}^n \mu(i)$$​

$1 \leq n \lt 2^{31}$

对于 $\varphi$ ，我们有 $\varphi(n) * 1 = id$，这里 $1$ 就是我们的 $g$，$id$ 就是我们的 $f*g$，这两个的求和是非常简单的。

对于 $\mu$，我们有 $\mu(n) * 1 = \epsilon(n) = [n = 1]$，这同样非常容易求和。

参考代码：

const int N = 5e6 + 10;
const int M = 5e3 + 10;

ll phi[N],maxn,maxm,n,mu[N],pri[N],tot,prephi[M],premu[M];
bool vis[N];

void init() {
    phi[1] = mu[1] = 1;
    for (int i = 2;i <= N-10;++i) {
        if (!vis[i]) pri[++tot] = i,mu[i] = -1,phi[i] = i - 1;
        for (int j = 1;j <= tot && pri[j] * i <= N - 10;++j) {
            vis[pri[j] * i] = 1;
            if (i % pri[j] == 0) {
                phi[i * pri[j]] = phi[i] * pri[j];
                break;
            }
            phi[i * pri[j]] = phi[i] * phi[pri[j]];
            mu[i * pri[j]] = -mu[i];
        }
    }
    for (int i = 1;i <= N-10;++i) mu[i] += mu[i-1],phi[i] += phi[i-1];
}//在O(n^c)的时间内进行预处理。

ll getphi(ll x) {
    if (x < N-10) return phi[x];
    ll now = n / x;
    if (prephi[now]) return prephi[now];
    //注意这边，我们记忆化的实际上是n / x的值，因为我们递归调用到的总是一个形如n/x的形式，所以这里的数组可以开的相对小一点
    ll ans = x * (x + 1) / 2;//\sum f*g
    for (ll l = 2,r;l <= x;l = r + 1) {
        r = x / (x / l);
        ans -= (r - l + 1) * getphi(x / l);//-\sum g
    }
    return prephi[now] = ans;//记得记忆化
}

ll getmu(ll x) {
    if (x < N-10) return mu[x];
    ll now = n / x;
    if (premu[now]) return premu[now];
    ll ans = 1;
    for (ll l = 2,r;l <= x;l = r + 1) {
        r = x / (x / l);
        ans -= (r - l + 1) * getmu(x / l);
    }
    return premu[now] = ans;
}

P3768 简单的数学题

求

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

$n \leq 10^{10}$

$$
\begin{aligned}
&\sum_{i=1}^n\sum_{j=1}^n ij \gcd(i,j)\
&= \sum_{i=1}^n\sum_{j=1}^n\sum_{d=1}^nd ij [\gcd(i,j) = d]\
&= \sum_{i=1}^{n/d}\sum_{j=1}^{n/d}\sum_{d=1}^nd^3 ij [\gcd(i,j) = 1]\
&= \sum_{d = 1} ^ n d^3 \sum\limits_{i=1}^{n/d}\sum\limits_{j=1}^{n/d}ij \sum\limits_{k | i,k | j} \mu(k)\
&= \sum_{d = 1} ^ n d^3 \sum\limits_{k} k^2\mu(k) \sum\limits_{i=1}^{n/dk}\sum\limits_{j=1}^{n/dk}ij \
&= \sum_{d = 1} ^ n d^3 \sum\limits_{k} k^2\mu(k) (\sum_{i=1}^{n/dk} i) ^2 \
&\text{Let sum(x)} = \sum\limits_{i=1}^x i,T = dk\
&= \sum_{T=1}^n sum(\frac{n}{T})^2 T^2\sum_{d|T} d\mu(\frac{T}{d})
\end{aligned}
$$

那么我们要求的就是 $T^2\sum_{d|T}d\mu(\frac{T}{d})$，然后我们发现后面这个东西再 Dirchlet 卷积一下有 $\sum_{d|T}d\mu(\frac{T}{d}) = \varphi(T)$，所以就变成了 $T^2\varphi(T)$ ，太神奇啦家人们！

那么我们总的式子就变成了：
$$
\sum_{T=1}^n sum(\frac{n}{T})^2 T^2\varphi(T)
$$
别急，按杜教筛套路继续推下去，我们令 $f(i) = i^2 \varphi(i)$ ， $F(i) = \sum\limits_{i=1}^n f(i)$

我们来考虑这个 $f$ 怎么卷比较可做！因为有 $\varphi * 1 = id$，所以我们直接设 $g(x) = x^2$ ，这样我们就易得 $f * g (x) = x ^3$。这就好起来了！所以我们的 $F$ 就好做了：
$$
\begin{aligned}
g(1)F(n) &= \sum\limits_{i = 1}^n f*g(i) -\sum\limits_{i=2}^n g(i)\times F(\lfloor \frac{n}{i}\rfloor) \
F(n) &= \sum\limits_{i=1}^n i^3 - \sum\limits_{i=2}^n i^2 \times F(\lfloor \frac{n}{i}\rfloor)
\end{aligned}
$$
这个是典型的杜教筛式，所以我们就数论分块里面再套一个杜教筛就OK了。

代码：

const int N = 8000000 + 10;
const int M = 1e5 + 10;

ll phi[N],maxn,maxm,n,mu[N],pri[N],tot,premu[M];
ll mod;
bool vis[N];
unordered_map <ll,ll> prephi;

void init() {
    phi[1] = mu[1] = 1;
    for (int i = 2;i <= N-10;++i) {
        if (!vis[i]) pri[++tot] = i,mu[i] = -1,phi[i] = i - 1;
        for (int j = 1;j <= tot && pri[j] * i <= N - 10;++j) {
            vis[pri[j] * i] = 1;
            if (i % pri[j] == 0) {
                phi[i * pri[j]] = phi[i] * pri[j];
                break;
            }
            phi[i * pri[j]] = phi[i] * phi[pri[j]];
            mu[i * pri[j]] = -mu[i];
        }
    }
    for (ll i = 1;i <= N-10;++i) phi[i] = (phi[i-1] + phi[i] * i % mod * i  % mod) % mod;
}


ll inv2,inv6;

ll sum2(ll x) {
    x %= mod;
    return x * (x + 1) % mod * inv2 % mod; 
}

ll sum3(ll x) {
    x %= mod;
    return x * (x + 1) % mod * (2 * x % mod + 1) %mod *inv6 %mod;
}

ll dusieve(ll x) {
    //printf("%lld ",x);
    if (x < N-10) return phi[x];
    ll now = n / x;
    if (prephi[x]) return prephi[x];
    ll ans = sum2(x) * sum2(x) % mod;
    for (ll l = 2,r;l <= x;l = r + 1) {
        r = x / (x / l);
        ans = ((ans - (((sum3(r) - sum3(l-1)) % mod + mod) % mod * dusieve(x / l))) % mod + mod) % mod;
    }
    return prephi[x] = (ans + mod) % mod;
}

ll ffpow(ll a,ll b) {
    a %= mod;
    ll ans = 1;
    for (;b;b >>= 1) {
        if (b & 1) ans = ans * a % mod;
        a = a * a % mod;
    }
    return ans;
}


signed main() {
    read(mod,n);
    inv2 = ffpow(2,mod - 2);inv6 = ffpow(6,mod - 2);
    init();
	ll ans = 0;
    for (ll l = 1,r;l <= n;l = r + 1) {
        r = n / (n / l);
        ll a1 = sum2(n / l) * sum2(n / l) % mod;
        ll a2 = ((dusieve(r) - dusieve(l - 1)) % mod + mod) % mod;
        ans = (ans + (a1 * a2) % mod) % mod;
    }
    if (ans < 0) ans += mod;
    printf("%lld\n",ans % mod);
    return 0;
}