杜教筛

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

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

杜教筛

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

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

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

$$egin{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} ight floor) \
&= sum _{d=2}^n g(d)S(left lfloor frac{n}{d} ight floor) + 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} ight floor)$.

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

莫比乌斯函数前缀和

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

$$egin{aligned}
herefore S_1(n) &= sum_{i=1}^n varepsilon (i) - sum_{i=2}^nS_1(left lfloor frac{n}{i} ight floor) \
&= 1 -  sum_{i=2}^nS_1(left lfloor frac{n}{i} ight floor)
end{aligned}$$

观察到 $left lfloor frac{n}{i} ight floor$ 最多只有 $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 egin{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)leftlfloorfrac{n}{d} ight floor^{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}})$

使用杜教筛求解

$ecause varphi * 1=I D$

$ herefore  S(n)=frac{1}{2} n(n+1)-sum_{i=2}^{n} Sleft(leftlfloorfrac{n}{i} ight floor ight)$

代码实现

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

#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
", S_phi(n), S_mu(n));
  }
  return 0;
}
原文地址:https://www.cnblogs.com/lfri/p/11691349.html