zoukankan      html  css  js  c++  java
  • 杜教筛

    转载自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;
    }
  • 相关阅读:
    Azure ARM (2) 概览
    Azure SQL Database (22) 迁移部分数据到Azure Stretch Database
    Azure SQL Database (21) 将整张表都迁移到Azure Stretch Database里
    Azure SQL Database (20) 使用SQL Server 2016 Upgrade Advisor
    开源日志系统比较:scribe、chukwa、kafka、flume
    大型互联网架构概述
    可扩展Web架构与分布式系统
    linux TOP命令各参数详解【转载】
    MongoDB社区版本和企业版本差别
    基于Geoserver的业务分析(笔记)
  • 原文地址:https://www.cnblogs.com/lfri/p/11691349.html
Copyright © 2011-2022 走看看