杜教筛

转载自https://oi-wiki.org/math/du/

在莫比乌斯反演的题目中,往往要求出一些数论函数的前缀和,利用 杜教筛 可以快速求出这些前缀和。

杜教筛

求 $\displaystyle S(n)=\sum_{i=1}^n f(i)$

我们要想办法构造一个 $S(n)$ 关于 $S(\left \lfloor \frac{n}{i} \right \rfloor)$.

构造两个积性函数 $g$,使得 $h=f*g$ 的前缀和很好求

$$\begin{aligned}
\sum _{i=1}^n h(i) &= \sum_{i=1}^n \sum _{d|n} g(d)f(\frac{n}{d}) \\
&= \sum_{d=1}^n g(d)\cdot \sum_{i=1}^{\frac{n}{d}} f(i)   \\
&= \sum_{d=1}^n g(d)S(\left \lfloor \frac{n}{d} \right \rfloor) \\
&= \sum _{d=2}^n g(d)S(\left \lfloor \frac{n}{d} \right \rfloor) + g(1)S(n)
\end{aligned}$$

因此 $g(1)S(n) = \sum_{i=1}^nh(i) - \sum _{d=2}^n g(d)S(\left \lfloor \frac{n}{d} \right \rfloor)$.

那么假设我们可以快速对 $\sum_{i=1}^n(f*g)(i)$ 求和,并用数论分块求解 $\sum_{i=2}^n g(i)S(\frac{n}{i})$ 就可以在较短时间内求得 $g(1)S(n)$.

莫比乌斯函数前缀和

前面证明过如下卷积式:$\mu * 1 = \varepsilon$

$$\begin{aligned}
\therefore S_1(n) &= \sum_{i=1}^n \varepsilon (i) - \sum_{i=2}^nS_1(\left \lfloor \frac{n}{i} \right \rfloor) \\
&= 1 -  \sum_{i=2}^nS_1(\left \lfloor \frac{n}{i} \right \rfloor)
\end{aligned}$$

观察到 $\left \lfloor \frac{n}{i} \right \rfloor$ 最多只有 $O(\sqrt n)$ 种取值,我们就可以应用整除分块(或称数论分块)来计算每一项的值了。

直接计算的时间复杂度为 $O(n^{\frac{3}{4}})$。考虑先线性筛预处理出前 $n^{\frac{2}{3}}$ 项,剩余部分的时间复杂度为 $O(\int _0^{n^\frac{1}{3}} \sqrt{\frac{n}{x}}dx) = O(n^{\frac{2}{3}})$.

对于较大值,需要用 map (或者unordered_map)存在其对应的值,方便以后使用时直接使用之前计算的结果。

欧拉函数前缀和

当然也可以用杜教筛求出 $\varphi (x)$ 的前缀和,但是更好的方法是应用莫比乌斯反演:

$\displaystyle \begin{array}
{l}{\sum_{i=1}^{n} \sum_{j=1}^{n} 1[\operatorname{gcd}(i, j)=1]=\sum_{i=1}^{n} \sum_{j=1}^{n} \sum_{d|i, d| j} \mu(d)} \\
{=\sum_{d=1}^{n} \mu(d)\left\lfloor\frac{n}{d}\right\rfloor^{2}}
\end{array}$

由于题目所求的是 $\sum_{i=1}^{n} \sum_{j=1}^{i} 1[\operatorname{gcd}(i, j)=1]$,所以我们排除掉 $i=j=1$ 的情况,并将结果除2即可。

观察到,只需求出莫比乌斯函数的前缀和,就可以快速计算出欧拉函数的前缀和了。时间复杂度为 $O(n^{\frac{2}{3}})$

使用杜教筛求解

$\because \varphi * 1=I D$

$\therefore  S(n)=\frac{1}{2} n(n+1)-\sum_{i=2}^{n} S\left(\left\lfloor\frac{n}{i}\right\rfloor\right)$

代码实现

//求莫比乌斯函数和欧拉函数的前缀和

#include <algorithm>
#include <cstdio>
#include <cstring>
#include <map>
using namespace std;
const int maxn = 2000010;
typedef long long ll;
ll T, n, pri[maxn], tot, mu[maxn], sum_mu[maxn];
bool vis[maxn];
map<ll, ll> mp_mu;  //可换成unordered_map
ll S_mu(ll x) {
  if (x < maxn) return sum_mu[x];
  if (mp_mu[x]) return mp_mu[x];
  ll ret = 1ll;
  for (ll i = 2, j; i <= x; i = j + 1) {
    j = x / (x / i);
    ret -= S_mu(x / i) * (j - i + 1);
  }
  return mp_mu[x] = ret;
}
ll S_phi(ll x) {
  ll ret = 0ll;
  for (ll i = 1, j; i <= x; i = j + 1) {
    j = x / (x / i);
    ret += (S_mu(j) - S_mu(i - 1)) * (x / i) * (x / i);
  }
  return ((ret - 1) >> 1) + 1;
}

void initMu()
{
  mu[1] = 1;
  for (int i = 2; i < maxn; i++) {
    if (!vis[i]) pri[++tot] = i, mu[i] = -1;
    for (int j = 1; j <= tot && i * pri[j] < maxn; j++) {
      vis[i * pri[j]] = true;
      if (i % pri[j] == 0)
      {
          mu[i * pri[j]] = 0;
          break;
      }
      else
      {
          mu[i * pri[j]] = -mu[i];
      }
    }
  }
  for (int i = 1; i < maxn; i++) sum_mu[i] = sum_mu[i - 1] + mu[i];
}

int main() {
  initMu();
  scanf("%lld", &T);
  while (T--) {
    scanf("%lld", &n);
    printf("%lld %lld\n", S_phi(n), S_mu(n));
  }
  return 0;
}

 

posted @ 2019-10-17 12:40  Rogn  阅读(241)  评论(0编辑  收藏  举报