杜教筛学习总结

前言

一直听说杜教筛非常nb，但关于它的学习一直被鸽= =但最近遇到太多数学题辣，所以不得不把这个坑填上了。

另外杜教筛可能需要一些前置知识，之前写过一篇关于莫比乌斯函数的，就顺便贴上来吧：传送门

正文

数论函数：我们平时遇到的一些特殊函数比如\(\varphi,\mu\)这种都属于数论函数。

积性函数：
定义：如果一个已知函数为数论函数，且\(f(1)=1\)，并且满足对于任意两个互质的\(p,q\)，有：\(f(p\cdot q)=f(p)\cdot f(q)\)，那么则称\(f\)为积性函数。特殊地，若对于任意两个数\(p,q\)，都有\(f(p\cdot q)=f(p)\cdot f(q)\)，那么就称这类函数为完全积性函数。

常见积性函数：

\(\mu(n)\)——莫比乌斯函数，这是十分常见，需要对其性质有一定了解；
\(\varphi(n)\)——欧拉函数，定义为\(1\)~\(n\)中与\(n\)互质的数的个数；
\(d(n)\)——约数个数，表示约数的个数。
(完全积性函数):
\(\varepsilon(n)\)——元函数，当\(n=1\)时其值为1，其余时候其值为0；
\(I(n)\)——恒等函数，不论\(n\)为何值，其值恒等于1；
\(id(n)\)——单位函数，\(id(n)=n\)；

这几个完全积性函数看似较为简单，其实有着很大的用处。它们完全积性的性质也是很好证明的。

狄利克雷卷积：
定义：两个函数\(f(n),g(n)\)的狄利克雷卷积记为\((f*g)(n)\)，表达式等价于\(\sum_{d|n}f(d)g(\frac{n}{d})\)。前面的括号代表卷积对象，后面的括号代表卷积的范围。
性质：
狄利克雷卷积运算具有一些很有用的性质：

交换率：\(f*g=g*f\);
结合律：\((f*g)*h=f*(g*h)\);
分配律：\(f*(g+h)=f*g+f*h\).

另外还有比较重要的一个性质：
若\(f,g\)都为积性函数，那么\(f*g\)也为积性函数。
在杜教筛中，通过卷积来判断积性函数，是很常用的一个技巧。
另外对于我们上面所说的元函数\(\varepsilon\)，有\(\varepsilon*f=f\)，即元函数类似于单位1在乘法中的作用。

应用：

对于莫比乌斯函数\(\mu\)来说，有\(\mu*I=\varepsilon\)；
对于欧拉函数\(\varphi\)来说，有\(\varphi*I=id\).

以上两个表达式也就是这两个函数最常用的性质之一。通过它们也可以来证明一些东西：
比如对于莫比乌斯反演形式一来说，已知\(F=f*I\)，我们将其两边同时卷上\(\mu\)，那么等式就为\(F*mu=f*\varepsilon=f\)，即\(f=F*\mu\)。
是不是感觉证明一些就容易许多了？
还有对于\(\varphi*I=id\)，两边同时卷上\(\mu\)，那么就有\(\varphi=id*\mu=\sum_{d|n}d\mu(\frac{n}{d})\)。这个式子也是很有用的。
总之，熟练掌握\(\mu,\varphi\)以及狄利克雷卷积的性质是很有好处的。

杜教筛

前面都是一些前置知识，接下来就步入正题——杜教筛啦。
杜教筛主要就是用来在压线性时间\(O(n^{\frac{2}{3}})\)内求出积性函数的前缀和。是不是感觉很厉害？好吧，其实杜教筛也并非那么的神奇，往下看就知道了~

现在对于积性函数\(f\)，要求\(\sum_{i=1}^nf(i)\)，我们将其记为\(S(n)\)。
构造\(h=f*g\)，那么有：

\[\begin{aligned} \sum_{i=1}^nh(i)=&\sum_{i=1}^nf*g\\ =&\sum_{i=1}^n\sum_{d|i}f(\frac{i}{d})g(d)\\ =&\sum_dg(d)\sum_{d|i}f(\frac{i}{d})\\ =&\sum_dg(d)S(\lfloor\frac{n}{d}\rfloor)\\ \end{aligned} \]

现在把第一项单独提出来，那么就有：

\[\sum_{i=1}^nh(i)=g(1)\cdot S(n)+\sum_{d=2}^ng(d)S(\lfloor\frac{n}{d}\rfloor) \]

式子等价于：

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

如果我们能够找到一个积性函数\(g\)，使得\(\sum_{i=1}^n(f*g)(i)\)这部分的前缀和能够快速求，而后面可以直接整除分块，那就能够起到加速的作用的。

举几个常见例子来说明吧：

求\(\sum_{i=1}^n\mu(i)\)。

我们知道:\(\mu*I=\varepsilon\)，而\(\varepsilon\)这个函数足够简单！是能够快速求和的，所以我们令\(g=\varepsilon\)，那么求和式子就能化为：

\[S(n)=1-\sum_{d=2}^nS(\lfloor\frac{n}{d}\rfloor) \]

对于后部分，整除分块处理即可。

求\(\sum_{i=1}^n\varphi(i)\)。

还是利用\(\varphi\)的性质，我们知道当它卷上\(I\)时，能够变成一个较简单的函数：\(id\)。这个是支持快速求和的，那么我们直接取\(g=I\)即可。

求\(\sum_{i=1}^ni\varphi(i)\)。

这个貌似就没有那么好处理了，假设现在有个\(g\)，我们将其卷积的形式表示出来：\(\sum_{d|n}(d\varphi(d))g(\frac{n}{d})\)。
注意观察这个式子，当\(g=id\)时，多处来的一个\(d\)能够被消去，并且剩余一个\(n\)可以提到外面去。
那么自然就会想到取\(g=id\)啦，这时\((f*g)(i)=i\sum_{d|n}\varphi(d)=i\)。
是不是感觉挺奇妙的，式子也变的很简单了。

求\(\sum_{i=1}^ni^2\varphi(i)\)。

跟上面同样的思路，先把卷积形式写出来，有：\(\sum_{d|n}d^2\varphi(d)g(\frac{n}{d})\)。
同样的，想到令\(g=id\)，就有：\(\sum_{i=1}^nf*g=\sum_{i=1}^ni^2\)，其实就是和上面同样的套路~

反正对于寻找\(g\)这一步，需要一定的观察力以及耐心，实在不行一个一个试就行了。

具体实现的话可以递归来实现，\(getsum(n)\)能够得到\(S(n)\)，然后递归处理即可。
证明的话可以用主定理来证，直接上杜教筛的话复杂度是\(O( n^{\frac{3}{4}})\)的，但我们可以先预处理出\(n^{\frac{2}{3}}\)的前缀和，那样复杂度就是\(O(n^{\frac{2}{3}})\)了。
可以加上一个记忆化的操作，可以更快一点，尤其是对于多组数据来说。
下面附上一个模板题的代码：
洛谷P4213

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = 5e6 + 5;
int mu[N], p[N];
ll phi[N];
bool chk[N];
unordered_map <int, ll> mp1, mp2;
void init() {
    mu[1] = phi[1] = 1;
    int cnt = 0, k = N - 1;
    for(int i = 2; i <= k; i++) {
        if(!chk[i]) p[++cnt] = i, mu[i] = -1, phi[i] = i - 1;
        for(int j = 1; j <= cnt && i * p[j] <= k; j++) {
            chk[i * p[j]] = 1;
            if(i % p[j] == 0) {mu[i * p[j]] = 0; phi[i * p[j]] = phi[i] * p[j]; break;}
            mu[i * p[j]] = -mu[i]; phi[i * p[j]] = phi[i] * (p[j] - 1);
        }
    }
    for(int i = 1; i < N; i++) mu[i] += mu[i - 1], phi[i] += phi[i - 1];
}
ll djs_mu(int n) {
    if(n <= 5000000) return mu[n];
    if(mp1[n]) return mp1[n];
    ll ans = 1;
    for(int i = 2, j; i <= n; i = j + 1) {
        j = n / (n / i);
        ans -= (j - i + 1) * djs_mu(n / i);
    }
    return mp1[n] = ans;
}
ll djs_phi(int n) {
    if(n <= 5000000) return phi[n];
    if(mp2[n]) return mp2[n];
    ll ans = 1ll * (n + 1) * n / 2;
    for(int i = 2, j; i <= n; i = j + 1) {
        j = n / (n / i);
        ans -= (j - i + 1) * djs_phi(n / i);
    }
    return mp2[n] = ans;
}
int n, T;
int main() {
    init();
    cin >> T;
    while(T--) {
        cin >> n;
        ll ans1 = djs_mu(n);
        ll ans2 = djs_phi(n);
        cout << ans2 << ' ' << ans1 << '\n';
    }
    return 0;
}

再举个例子吧~

HDU6706 huntian oy
题意：
定义\(f(n,a,b)=\sum_{i=1}^n\sum_{j=1}^igcd(i^a-j^a,i^b-j^b)[gcd(i,j)=1]\% 10^9+7\)。
先有\(T\)组询问，每组询问给出\(n,a,b\)，求\(f(n,a,b)\)。

思路：
先有一个公式如下：
若\(a\geq b\)，且\(a,b\)互质，则有：

\[\begin{aligned} &gcd(a^n-b^n,a^m-b^m)\\ =&gcd(a^m-b^m,a^{n\%m}-b^{n\%m})=a^{gcd(n,m)}-b^{gcd(n,m)} \end{aligned} \]

我也不知道这咋证= =
反正有了这个式子之后就好推了：

\[\begin{aligned} f(n,a,b)=&\sum_{i=1}^n \sum_{j=1}^i i-j[gcd(i,j)=1]\\ =&\sum_{i=1}^n \sum_{j=1}^i i[gcd(i,j)=1]-\sum_{i=1}^n \sum_{j=1}^i j[gcd(i,j)=1]\\ =&\sum_{i=1}^ni\varphi(i)-\sum_{i=1}^n\frac{i\varphi(i)+[i=1]}{2} \end{aligned} \]

至于为什么后面那部分可以化成这样，这涉及到\(gcd(n,m)=gcd(n-m,n)\)，也就是说，与之互素的数是成对出现的，并且其和为\(n\)。当然\(n=1\)的情况除外。
之后式子就可以化为这样：

\[f(n,a,b)=\frac{1}{2}\sum_{i=1}^ni\varphi(i)-1 \]

之后直接上杜教筛就是了。

Code

#include <bits/stdc++.h>
#define heyuhhh ok
using namespace std;
typedef long long ll;
const int N = 5e6 + 5, MOD = 1e9 + 7, inv6 = 166666668;
int p[N];
ll phi[N];
bool chk[N];
unordered_map <int, int> mp;
inline void init() {
    phi[1] = 1;
    int cnt = 0, k = N - 1, i;
    for(i = 2; i <= k; ++i) {
        if(!chk[i]) p[++cnt] = i, phi[i] = i - 1;
        for(int j = 1; j <= cnt && i * p[j] <= k; j++) {
            chk[i * p[j]] = 1;
            if(i % p[j] == 0) {phi[i * p[j]] = phi[i] * p[j]; break;}
            phi[i * p[j]] = phi[i] * (p[j] - 1);
        }
    }
    for(i = 1; i < N; ++i) {
        phi[i] = phi[i - 1] + 1ll * i * phi[i];
        if(phi[i] >= MOD) phi[i] %= MOD;
    }
}

int djs_phi(int n) {
    if(n <= 5000000) return phi[n];
    if(mp[n]) return mp[n];
    int ans = 1ll * (n + 1) * n % MOD * (2ll * n + 1) % MOD * inv6 % MOD;
    for(register int i = 2, j; i <= n; i = j + 1) {
        j = n / (n / i);
        ans -= ((ll(j - i + 1) * (j + i)) >> 1) % MOD * (ll)djs_phi(n / i) % MOD;
        ans %= MOD;
    }
    if(ans < 0) ans += MOD;
    return mp[n] = ans;
}
int n, a, b, T;
int main() {
#ifdef heyuhhh
    freopen("input.in", "r", stdin);
#else
    ios::sync_with_stdio(false); cin.tie(0);
#endif
    init();
    cin >> T;
    while(T--) {
        cin >> n >> a >> b;
        int ans = djs_phi(n) - 1;
        ans = (1ll * ans * (MOD + 1) / 2) % MOD;
        cout << ans << '\n';
    }
    return 0;
}

总结

昨晚匆匆忙忙把这篇博客写完了=，=，感觉还是有很多东西没有写好，emmm以后再慢慢改吧。

posted @ 2019-08-25 21:43 heyuhhh 阅读(245) 评论(0) 编辑收藏举报

刷新页面返回顶部

heyuhhh

劳谦君子，有终吉。

杜教筛学习总结

杜教筛学习总结

前言

正文

杜教筛

总结

公告