跳转至

Min_25 篩

定義

從此種篩法的思想方法來説，其又被稱為「Extended Eratosthenes Sieve」。

由於其由 Min_25 發明並最早開始使用，故稱「Min_25 篩」。

性質

其可以在 \(O\left(\frac{n^{\frac{3}{4}}}{\log{n}}\right)\) 或 \(\Theta\left(n^{1 - \epsilon}\right)\) 的時間複雜度下解決一類 積性函數 的前綴和問題。

要求：\(f(p)\) 是一個關於 \(p\) 的項數較少的多項式或可以快速求值；\(f(p^{c})\) 可以快速求值。

記號

  • 如無特別説明，本節中所有記為 \(p\) 的變量的取值集合均為全體質數。
  • \(x / y := \left\lfloor\frac{x}{y}\right\rfloor\)
  • \(\operatorname{isprime}(n) := [ |\{d : d \mid n\}| = 2 ]\)，即 \(n\) 為質數時其值為 \(1\)，否則為 \(0\)。
  • \(p_{k}\)：全體質數中第 \(k\) 小的質數（如：\(p_{1} = 2, p_{2} = 3\)）。特別地，令 \(p_{0} = 1\)。
  • \(\operatorname{lpf}(n) := [1 < n] \min\{p : p \mid n\} + [1 = n]\)，即 \(n\) 的最小質因數。特別地，\(n=1\) 時，其值為 \(1\)。
  • \(F_{\mathrm{prime}}(n) := \sum_{2 \le p \le n} f(p)\)
  • \(F_{k}(n) := \sum_{i = 2}^{n} [p_{k} \le \operatorname{lpf}(i)] f(i)\)

解釋

觀察 \(F_{k}(n)\) 的定義，可以發現答案即為 \(F_{1}(n) + f(1) = F_{1}(n) + 1\)。

考慮如何求出 \(F_{k}(n)\)。通過枚舉每個 \(i\) 的最小質因子及其次數可以得到遞推式：

\[ \begin{aligned} F_{k}(n) &= \sum_{i = 2}^{n} [p_{k} \le \operatorname{lpf}(i)] f(i) \\ &= \sum_{\substack{k \le i \\ p_{i}^{2} \le n}} \sum_{\substack{c \ge 1 \\ p_{i}^{c} \le n}} f\left(p_{i}^{c}\right) ([c > 1] + F_{i + 1}\left(n / p_{i}^{c}\right)) + \sum_{\substack{k \le i \\ p_{i} \le n}} f(p_{i}) \\ &= \sum_{\substack{k \le i \\ p_{i}^{2} \le n}} \sum_{\substack{c \ge 1 \\ p_{i}^{c} \le n}} f\left(p_{i}^{c}\right) ([c > 1] + F_{i + 1}\left(n / p_{i}^{c}\right)) + F_{\mathrm{prime}}(n) - F_{\mathrm{prime}}(p_{k - 1}) \\ &= \sum_{\substack{k \le i \\ p_{i}^{2} \le n}} \sum_{\substack{c \ge 1 \\ p_{i}^{c + 1} \le n}} \left(f\left(p_{i}^{c}\right) F_{i + 1}\left(n / p_{i}^{c}\right) + f\left(p_{i}^{c + 1}\right)\right) + F_{\mathrm{prime}}(n) - F_{\mathrm{prime}}(p_{k - 1}) \end{aligned} \]

最後一步推導基於這樣一個事實：對於滿足 \(p_{i}^{c} \le n < p_{i}^{c + 1}\) 的 \(c\)，有 \(p_{i}^{c + 1} > n \iff n / p_{i}^{c} < p_{i} < p_{i + 1}\)，故 \(F_{i + 1}\left(n / p_{i}^{c}\right) = 0\)。
其邊界值即為 \(F_{k}(n) = 0 (p_{k} > n)\)。

假設現在已經求出了所有的 \(F_{\mathrm{prime}}(n)\)，那麼有兩種方式可以求出所有的 \(F_{k}(n)\)：

  1. 直接按照遞推式計算。
  2. 從大到小枚舉 \(p\) 轉移，僅當 \(p^{2} < n\) 時轉移增加值不為零，故按照遞推式後綴和優化即可。

現在考慮如何計算 \(F_{\mathrm{prime}}{(n)}\)。
觀察求 \(F_{k}(n)\) 的過程，容易發現 \(F_{\mathrm{prime}}\) 有且僅有 \(1, 2, \dots, \left\lfloor\sqrt{n}\right\rfloor, n / \sqrt{n}, \dots, n / 2, n\) 這 \(O(\sqrt{n})\) 處的點值是有用的。
一般情況下，\(f(p)\) 是一個關於 \(p\) 的低次多項式，可以表示為 \(f(p) = \sum a_{i} p^{c_{i}}\)。
那麼對於每個 \(p^{c_{i}}\)，其對 \(F_{\mathrm{prime}}(n)\) 的貢獻即為 \(a_{i} \sum_{2 \le p \le n} p^{c_{i}}\)。
分開考慮每個 \(p^{c_{i}}\) 的貢獻，問題就轉變為了：給定 \(n, s, g(p) = p^{s}\)，對所有的 \(m = n / i\)，求 \(\sum_{p \le m} g(p)\)。

Notice：\(g(p) = p^{s}\) 是完全積性函數！

於是設 \(G_{k}(n) := \sum_{i = 1}^{n} \left[p_{k} < \operatorname{lpf}(i) \lor \operatorname{isprime}(i)\right] g(i)\)，即埃篩第 \(k\) 輪篩完後剩下的數的 \(g\) 值之和。
對於一個合數 \(x\)，必定有 \(\operatorname{lpf}(x) \le \sqrt{x}\)，則 \(F_{\mathrm{prime}} = G_{\left\lfloor\sqrt{n}\right\rfloor}\)，故只需篩到 \(G_{\left\lfloor\sqrt{n}\right\rfloor}\) 即可。
考慮 \(G\) 的邊界值，顯然為 \(G_{0}(n) = \sum_{i = 2}^{n} g(i)\)。（還記得嗎？特別約定了 \(p_{0} = 1\)）
對於轉移，考慮埃篩的過程，分開討論每部分的貢獻，有：

  1. 對於 \(n < p_{k}^{2}\) 的部分，\(G\) 值不變，即 \(G_{k}(n) = G_{k - 1}(n)\)。
  2. 對於 \(p_{k}^{2} \le n\) 的部分，被篩掉的數必有質因子 \(p_{k}\)，即 \(-g(p_{k}) G_{k - 1}(n / p_{k})\)。
  3. 對於第二部分，由於 \(p_{k}^{2} \le n \iff p_{k} \le n / p_{k}\)，故會有 \(\operatorname{lpf}(i) < p_{k}\) 的 \(i\) 被減去。這部分應當加回來，即 \(g(p_{k}) G_{k - 1}(p_{k - 1})\)。

則有：

\[ G_{k}(n) = G_{k - 1}(n) - \left[p_{k}^{2} \le n\right] g(p_{k}) (G_{k - 1}(n / p_{k}) - G_{k - 1}(p_{k - 1})) \]

複雜度分析

對於 \(F_{k}(n)\) 的計算，其第一種方法的時間複雜度被證明為 \(O\left(n^{1 - \epsilon}\right)\)（見 zzt 集訓隊論文 2.3）；
對於第二種方法，其本質即為洲閣篩的第二部分，在洲閣論文中也有提及（6.5.4），其時間複雜度被證明為 \(O\left(\frac{n^{\frac{3}{4}}}{\log{n}}\right)\)。

對於 \(F_{\mathrm{prime}}(n)\) 的計算，事實上，其實現與洲閣篩第一部分是相同的。
考慮對於每個 \(m = n / i\)，只有在枚舉滿足 \(p_{k}^{2} \le m\) 的 \(p_{k}\) 轉移時會對時間複雜度產生貢獻，則時間複雜度可估計為：

\[ \begin{aligned} T(n) &= \sum_{i^{2} \le n} O\left(\pi\left(\sqrt{i}\right)\right) + \sum_{i^{2} \le n} O\left(\pi\left(\sqrt{\frac{n}{i}}\right)\right) \\ &= \sum_{i^{2} \le n} O\left(\frac{\sqrt{i}}{\ln{\sqrt{i}}}\right) + \sum_{i^{2} \le n} O\left(\frac{\sqrt{\frac{n}{i}}}{\ln{\sqrt{\frac{n}{i}}}}\right) \\ &= O\left(\int_{1}^{\sqrt{n}} \frac{\sqrt{\frac{n}{x}}}{\log{\sqrt{\frac{n}{x}}}} \mathrm{d} x\right) \\ &= O\left(\frac{n^{\frac{3}{4}}}{\log{n}}\right) \end{aligned} \]

對於空間複雜度，可以發現不論是 \(F_{k}\) 還是 \(F_{\mathrm{prime}}\)，其均只在 \(n / i\) 處取有效點值，共 \(O(\sqrt{n})\) 個，僅記錄有效值即可將空間複雜度優化至 \(O(\sqrt{n})\)。

首先，通過一次數論分塊可以得到所有的有效值，用一個大小為 \(O(\sqrt{n})\) 的數組 \(\text{lis}\) 記錄。對於有效值 \(v\)，記 \(\text{id}(v)\) 為 \(v\) 在 \(\text{lis}\) 中的下標，易得：對於所有有效值 \(v\)，\(\text{id}(v) \le \sqrt{n}\)。

然後分開考慮小於等於 \(\sqrt{n}\) 的有效值和大於 \(\sqrt{n}\) 的有效值：對於小於等於 \(\sqrt{n}\) 的有效值 \(v\)，用一個數組 \(\text{le}\) 記錄其 \(\text{id}(v)\)，即 \(\text{le}_v = \text{id}(v)\)；對於大於 \(\sqrt{n}\) 的有效值 \(v\)，用一個數組 \(\text{ge}\) 記錄 \(\text{id}(v)\)，由於 \(v\) 過大所以藉助 \(v' = n / v < \sqrt{n}\) 記錄 \(\text{id}(v)\)，即 \(\text{ge}_{v'} = \text{id}(v)\)。

這樣，就可以使用兩個大小為 \(O(\sqrt{n})\) 的數組記錄所有有效值的 \(\text{id}\) 並 \(O(1)\) 查詢。在計算 \(F_{k}\) 或 \(F_{\mathrm{prime}}\) 時，使用有效值的 \(\text{id}\) 代替有效值作為下標，即可將空間複雜度優化至 \(O(\sqrt{n})\)。

過程

對於 \(F_{k}(n)\) 的計算，我們實現時一般選擇實現難度較低的第一種方法，其在數據規模較小時往往比第二種方法的表現要好；

對於 \(F_{\mathrm{prime}}(n)\) 的計算，直接按遞推式實現即可。

對於 \(p_{k}^{2} \le n\)，可以用線性篩預處理出 \(s_{k} := F_{\mathrm{prime}}(p_{k})\) 來替代 \(F_{k}\) 遞推式中的 \(F_{\mathrm{prime}}(p_{k - 1})\)。
相應地，\(G\) 遞推式中的 \(G_{k - 1}(p_{k - 1}) = \sum_{i = 1}^{k - 1} g(p_{i})\) 也可以用此方法預處理。

用 Extended Eratosthenes Sieve 求 積性函數 \(f\) 的前綴和時，應當明確以下幾點：

  • 如何快速（一般是線性時間複雜度）篩出前 \(\sqrt{n}\) 個 \(f\) 值；
  • \(f(p)\) 的多項式表示；
  • 如何快速求出 \(f(p^{c})\)。

明確上述幾點之後按順序實現以下幾部分即可：

  1. 篩出 \([1, \sqrt{n}]\) 內的質數與前 \(\sqrt{n}\) 個 \(f\) 值；
  2. 對 \(f(p)\) 多項式表示中的每一項篩出對應的 \(G\)，合併得到 \(F_{\mathrm{prime}}\) 的所有 \(O(\sqrt{n})\) 個有用點值；
  3. 按照 \(F_{k}\) 的遞推式實現遞歸，求出 \(F_{1}(n)\)。

例題

求莫比烏斯函數的前綴和

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

易知 \(f(p) = -1\)。則 \(g(p) = -1, G_{0}(n) = \sum_{i = 2}^{n} g(i) = -n + 1\)。
直接篩即可得到 \(F_{\mathrm{prime}}\) 的所有 \(O(\sqrt{n})\) 個所需點值。

求歐拉函數的前綴和

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

首先易知 \(f(p) = p - 1\)。
對於 \(f(p)\) 的一次項 \((p)\)，有 \(g(p) = p, G_{0}(n) = \sum_{i = 2}^{n} g(i) = \frac{(n + 2) (n - 1)}{2}\)；
對於 \(f(p)\) 的常數項 \((-1)\)，有 \(g(p) = -1, G_{0}(n) = \sum_{i = 2}^{n} g(i) = -n + 1\)。
篩兩次加起來即可得到 \(F_{\mathrm{prime}}\) 的所有 \(O(\sqrt{n})\) 個所需點值。

「LOJ #6053」簡單的函數

給定 \(f(n)\)：

\[ f(n) = \begin{cases} 1 & n = 1 \\ p \operatorname{xor} c & n = p^{c} \\ f(a)f(b) & n = ab \land a \perp b \end{cases} \]

易知 \(f(p) = p - 1 + 2[p = 2]\)。則按照篩 \(\varphi\) 的方法篩，對 \(2\) 討論一下即可。
此處給出一種 C++ 實現：

參考代碼 
  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
/* 「LOJ #6053」简单的函数 */
#include <algorithm>
#include <cmath>
#include <cstdio>

const int maxs = 200000;  // 2sqrt(n)
const int mod = 1000000007;

template <typename x_t, typename y_t>
void inc(x_t &x, const y_t &y) {
  x += y;
  (mod <= x) && (x -= mod);
}

template <typename x_t, typename y_t>
void dec(x_t &x, const y_t &y) {
  x -= y;
  (x < 0) && (x += mod);
}

template <typename x_t, typename y_t>
int sum(const x_t &x, const y_t &y) {
  return x + y < mod ? x + y : (x + y - mod);
}

template <typename x_t, typename y_t>
int sub(const x_t &x, const y_t &y) {
  return x < y ? x - y + mod : (x - y);
}

template <typename _Tp>
int div2(const _Tp &x) {
  return ((x & 1) ? x + mod : x) >> 1;
}

// 以上目的均为防负数和取模
template <typename _Tp>
long long sqrll(const _Tp &x) {  // 平方函数
  return (long long)x * x;
}

int pri[maxs / 7], lpf[maxs + 1], spri[maxs + 1], pcnt;

void sieve(const int &n) {
  for (int i = 2; i <= n; ++i) {
    if (lpf[i] == 0) {  // 记录质数
      lpf[i] = ++pcnt;
      pri[lpf[i]] = i;
      spri[pcnt] = sum(spri[pcnt - 1], i);  // 前缀和
    }
    for (int j = 1, v; j <= lpf[i] && (v = i * pri[j]) <= n; ++j) lpf[v] = j;
  }
}

long long global_n;
int lim;
int le[maxs + 1],  // x <= \sqrt{n}
    ge[maxs + 1];  // x > \sqrt{n}
#define idx(v) (v <= lim ? le[v] : ge[global_n / v])

int G[maxs + 1][2], Fprime[maxs + 1];
long long lis[maxs + 1];
int cnt;

void init(const long long &n) {
  for (long long i = 1, j, v; i <= n; i = n / j + 1) {
    j = n / i;
    v = j % mod;
    lis[++cnt] = j;
    (j <= lim ? le[j] : ge[global_n / j]) = cnt;
    G[cnt][0] = sub(v, 1ll);
    G[cnt][1] = div2((long long)(v + 2ll) * (v - 1ll) % mod);
  }
}

void calcFprime() {
  for (int k = 1; k <= pcnt; ++k) {
    const int p = pri[k];
    const long long sqrp = sqrll(p);
    for (int i = 1; lis[i] >= sqrp; ++i) {
      const long long v = lis[i] / p;
      const int id = idx(v);
      dec(G[i][0], sub(G[id][0], k - 1));
      dec(G[i][1], (long long)p * sub(G[id][1], spri[k - 1]) % mod);
    }
  }
  /* F_prime = G_1 - G_0 */
  for (int i = 1; i <= cnt; ++i) Fprime[i] = sub(G[i][1], G[i][0]);
}

int f_p(const int &p, const int &c) {
  /* f(p^{c}) = p xor c */
  return p xor c;
}

int F(const int &k, const long long &n) {
  if (n < pri[k] || n <= 1) return 0;
  const int id = idx(n);
  long long ans = Fprime[id] - (spri[k - 1] - (k - 1));
  if (k == 1) ans += 2;
  for (int i = k; i <= pcnt && sqrll(pri[i]) <= n; ++i) {
    long long pw = pri[i], pw2 = sqrll(pw);
    for (int c = 1; pw2 <= n; ++c, pw = pw2, pw2 *= pri[i])
      ans +=
          ((long long)f_p(pri[i], c) * F(i + 1, n / pw) + f_p(pri[i], c + 1)) %
          mod;
  }
  return ans % mod;
}

int main() {
  scanf("%lld", &global_n);
  lim = sqrt(global_n);  // 上限

  sieve(lim + 1000);  // 预处理
  init(global_n);
  calcFprime();
  printf("%lld\n", (F(1, global_n) + 1ll + mod) % mod);

  return 0;
}
本页面最近更新：，更新历史
发现错误？想一起完善？ 在 GitHub 上编辑此页！
本页面贡献者：Marcythm, Xeonacid
本页面的全部内容在 CC BY-SA 4.0 和 SATA 协议之条款下提供，附加条款亦可能应用