# 質數篩法
質數是個非常深奧的數字，在這個單元會討論如何判斷一個正整數是否為質數
## 複習: 什麼是質數?
設兩個數字 $p,q$ 乘起來會變成另一個數 $N$。那麼就可以說 $p,q$ 是 $N$ 的因數。顯然 $N$ 有可能會有很多組因數 $(p,q)$，舉例像是:
$$
\begin{split}
12&=1\cdot 12\\
&=2\cdot 6\\
&=3\cdot 4\\
\end{split}
$$
這樣我們會說 $12$ 的因數有: $\{1,2,3,4,6,12\}$
透過上面的定義不難觀察到以下性質:
- 任何正整數 $N$ 的因數肯定有 $1$ 與 $N$
- 除了 $1,N$ 以外的其他因數會被包含在區間 $[2, N-1]$
有些數字是很特殊的存在，他們的因數都只有自己本身與 $1$，像是 $7$ 的因數僅有 $\{1,7\}$。這種數字我們稱為「質數」，因此我們可以觀察到另外一點:
- 質數在區間 $[2,N-1]$ 不會有任何因數
## 判斷質數
### 暴力作法判斷質數
我們回到我們的問題，如何判斷一個正整數 $N$ 是質數? 透過上面觀察到的性質，不難想出最簡單且直觀的方法，就是遍歷整個區間 $[2, N-1]$，每個數字 $i$ 都除除看，找餘數為何 (在程式中有個很方便的運算子 "%" 可以算餘數)。如果餘數是 $0$ 代表 $N$ 可以整除 $i$，且 $i$ 是 $N$ 的因數
```cpp=
bool isPrime(int N) {
for (int i = 2; i <= N - 1; i++) { // 檢查整個區間 [2, N - 1]
if(N % i == 0) {
return false;
}
}
return true;
}
```
這樣所花的時間顯然為 $O(N-2)=O(N)$。此方法實在太差，或許我們可以考慮縮小區間
### $O(\sqrt N)$ 作法判斷質數
觀察一下等式 $N=p\cdot q$，當 $p$ 能整除 $N$，那 $q$ 當然也能整除 $N$，因此沒必要重複判斷 $q$。反之，如果某整數 $i$ 已經確認不是因數，則 $\tfrac{N}{i}$ 當然也不是因數。再者，因為判斷 $i$ 時已經判斷過 $\tfrac{N}{i}$，所以不需再判斷一次 $\tfrac{N}{i}$。根據等式 $N=p\cdot q$，隨著 $p$ 上升，相對應 $q$ 下降，兩者必然有交點 $p=q=\sqrt{N}$
根據以上，其實只需要判斷到 $[2,\sqrt{N}]$ 就好。我們可以假設要搜尋的整數為 $i$，令搜尋條件為 `i*i<=N`，然後來跑整個迴圈
```cpp=
bool isPrime(int N) {
for (int i = 2; i * i <= N; i++) {
if(N % i == 0) {
return false;
}
}
return true;
}
```
這裡不只用 `i<=sqrt(N)` 的原因是因為可能會有浮點數的問題，這種誤差常常會出問題。在寫程式的時候，浮點數能避開就避，因為很多數字光是儲存就一定會有誤差，更不用說計算時誤差會傳導
## 質數篩法 Sieve of Eratosthenes
這裡講的質數篩法是「埃拉托斯特尼篩法」(埃氏篩法)，也就是古希臘人用的那個方法。假設今天需要判斷的正整數總共有 $q$ 筆資料，那麼時間花費為 $O(q\sqrt{N})$。如果考慮預處理一個表格紀錄「每個數字是否為質數」那麼其實是可以把時間降到 $O(q)$
### 演算法
可以觀察到如果一個數字 $n$ 已經確定是質數，那麼對於所有正整數 $c$，$c\cdot n$ 都是合數，而剩下的就是質數
<center>
<img src="https://upload.wikimedia.org/wikipedia/commons/b/b9/Sieve_of_Eratosthenes_animation.gif" style="margin: 0 auto; display: block; width: 300px">
<br>
圖取自維基百科
</center>
在實作上，我們需要先建立一個陣列 `is_prime[]` 並開兩個迴圈:
- 外層找質數
- 內層將該質數的倍數都設為「非質數」
跟上面的方法一樣，外層其實可以只跑區間 $[2, N-1]$ 就好。而內層迴圈中，由於目前跑的數字的平方其實已經標記過了，所以從他開始繼續標記下去
```cpp=
int N;
vector<bool> is_prime(N + 1, true);
is_prime[0] = is_prime[1] = false;
void sieve() {
for (int i = 2; i * i <= N; i++) {
if (!is_prime[i]) continue;
for (int j = i * i; j <= N; j += i)
is_prime[j] = false;
}
}
```
### 複雜度分析
- 時間: $O(N\cdot \log \log \sqrt{N})$，推導跟調和級數有關，這邊就不證了
- 空間: $O(N)$
由於`is_prime[]` 是預處理的陣列，一旦建立好就不需要再建一次，所以當我們有 $q$ 比詢問時，時間會是 $O(q)$。要注意每次可以查詢的區間為 $[0, N]$，超過這個區間就沒辦法判斷了。有時候題目可能輸入的是更大的數字時 (例: $10^7$)，時間會爆掉，因此必須考慮其他方法
## 線性篩法
實際在競賽時用埃氏篩法就夠了，根據我的經驗如果再不行就直接砸 Miller-Rabin 演算法，在大多數時候都可以過。但因為我看過這個神奇的線性篩法 (又稱歐拉篩法) 所以還是再多提一點
### 演算法
之前我們的方法是先遍歷整個區間 $[2,N-1]$，找一個質數 $q$，然後篩掉整個區間內所有以 $q$ 為質因數的合數。但是這樣會造成一個問題，就是兩個不同質數可能會篩到同一個數字
我們可以先觀察一個性質:
- 每個合數 $x$ 都只有一組 $(i, p)$，使得 $x=i\cdot p$，且 $p$ 為 $x$ 最小質因數
根據這個性質，就可以得出以下演算法
```cpp=
vector<int> primes; // 紀錄質數有誰
vector<bool> is_prime(N + 5, true); // 每個數字是否為質數
void linear_sieve() {
is_prime[0] = is_prime[1] = false;
for (int i = 2; i <= N; i++) {
if (is_prime[i]) { // 還沒被篩掉就是質數
primes.push_back(i);
}
for (int p : primes) {
if (i * p > N) break; // 超過就 break
is_prime[i * p] = false; // 篩掉
if (i % p == 0) break; // 下面會講
}
}
}
```
### 正確性
在 `if (i % p == 0)` 以前的程式應該都很好懂，但是這個 `if` 判斷就不知道在幹嘛，這很正常。注意到第二層迴圈中，會遍歷目前找到的所有質數 `p`，若是 `break` 代表找到一個質數可以被整除，當然此質數也會是最小。具體來說這件事情可以分兩部分討論:
- 如果 `i % p == 0`: 這代表 $p$ 就是 $i$ 的最小質因數。$p$ 雖然是 $x = i \times p$ 的最小質因數，但下個質數肯定不是，所以可以 `break`
- 如果 `i % p != 0`: 說明 $p$ 比 $i$ 的任何質因數都小。因此，$p$ 絕對是 $x = i \times p$ 的最小質因數，可以放心標記
### 複雜度分析
因為這方法會跑過所有合數，且不會重複，因此是 $O(N)$
## 應用: 質因數分解
### 預處理: 線性篩
我們可以考慮將最小質因數記錄到 `lp[]` 陣列當中，剩下都一樣，原理剛剛說過了
```cpp=
vector<int> primes;
vector<bool> is_prime(N + 5, true);
vector<int> lp(N + 5, 0);
void linear_sieve() {
is_prime[0] = is_prime[1] = false;
for (int i = 2; i <= N; ++i) {
if (is_prime[i]) {
lp[i] = i;
primes.push_back(i);
}
for (int p : primes) {
if (i * p > N)
break;
is_prime[i * p] = false;
lp[i * p] = p;
if (i % p == 0)
break;
}
}
}
```
### 質因數分解
一個數除以自己當前最小的質因數，直到不能再除就停止。code 還蠻好懂得，就不多做說明
```cpp=
vector<int> factorize(int x) {
vector<int> fac;
while (x > 1) {
fac.push_back(lp[x]); // 每次提取最小質因數
x /= lp[x]; // 除掉它，繼續分解剩下的部分
}
return fac;
}
```
### 時間複雜度
**線性篩**: $O(N)$
**分解**: 每次 `while` 保守估計會除掉 $2$，因此大概 $O(\log x)$ 次會把 $x$ 降為 $1$ 然後結束迴圈
## 題目練習
[Zerojudge b513. 判斷質數-商競103](https://zerojudge.tw/ShowProblem?problemid=b513)
[Zerojudge a626. 6. Prime Directive](https://zerojudge.tw/ShowProblem?problemid=a626)
[Zerojudge d709. 判断质数（一）](https://zerojudge.tw/ShowProblem?problemid=d709)
[Zerojudge a007. 判斷質數](https://zerojudge.tw/ShowProblem?problemid=a007)
[Zerojudge d362. 10394 - Twin Primes](https://zerojudge.tw/ShowProblem?problemid=d362)
[Zerojudge a010. 因數分解](https://zerojudge.tw/ShowProblem?problemid=a010)
[Zerojudge a537. 10789 - Prime Frequency](https://zerojudge.tw/ShowProblem?problemid=a537)
[Zerojudge a740. 质因数之和](https://zerojudge.tw/ShowProblem?problemid=a740)
[AtCoder Beginner Contest 250D - 250-like Number](https://atcoder.jp/contests/abc250/tasks/abc250_d?lang=en)
[Codeforces 142B. T-primes](https://codeforces.com/problemset/problem/230/B)
[Zerojudge a121. 質數又來囉](https://zerojudge.tw/ShowProblem?problemid=a121)
[Zerojudge d120. 10699 - Count the factors](https://zerojudge.tw/ShowProblem?problemid=d120)
[CSES Common Divisors](https://cses.fi/problemset/task/1081)
[Zerojudge c668. 下一個質數](https://zerojudge.tw/ShowProblem?problemid=c668) (||哈哈||)
[Codeforces 560D. Almost All Divisors](https://codeforces.com/problemset/problem/1165/D)
[Zerojudge e655. 10852 - Less Prime](https://zerojudge.tw/ShowProblem?problemid=e655)
[Codeforces 900F. Vasilije Loves Number Theory](https://codeforces.com/problemset/problem/1878/F)
----
## 參考資料
- [CP Algorithms-Sieve of Eratosthenes](https://cp-algorithms.com/algebra/sieve-of-eratosthenes.html#linear-time-modification)
- [CP Algorithms-Linear Sieve](https://cp-algorithms.com/algebra/prime-sieve-linear.html)
- [OI Wiki-篩法](https://oi-wiki.org/math/number-theory/sieve/)
----
> [ShanC 程式競賽筆記](https://hackmd.io/@ShanC/B1ouGxqcC)
> 作者: ShanC
> 更新: 2026/1/2
Sign in with Wallet
Connect another wallet