2020 競技程設教材 HackMD Book

2019 Week 15: Number theory & Calculation

數大便是美

模反元素

數字們同餘

$n$ 後，對於加法、減法、乘法，其性質都與以往差不多，但對數字做除法卻不盡人意

注意，在同餘運算中只會有整數，有理數無理數等其他數不會出現

反元素指的是元素與其運算後為單位元素的元素
例如：

• 整數
  $a$ 的加法反元素為
  $-a$
• 整數
  $a$ 的乘法反元素為
  $a^{-1}$

而元素

$a$ 的模

$n$ 反元素為

$a^{-1}$ 滿足

$$a\cdot a^{-1}\equiv 1\mathrm{mod}n$$

根據 Bezout's Thm

• $ax+ny=1$ 在模
  $n$ 後有
  $ax\equiv 1\mathrm{mod}n$
• $ax+ny\ne 1$ 在模
  $n$ 後有
  $ax\not\equiv 1\mathrm{mod}n$

也就是說

$gcd(a,n)=1$ (互質)，

$a$ 才擁有模

$n$ 反元素

求出模反元素

• 當
  $n$ 為質數時：

根據費馬小定理有

$$a\cdot a^{n-2}\equiv 1\mathrm{mod}n$$
所以

$a^{n-2}$ 是

$a$ 的模反元素，可用快速冪 在

$O(\log n)$ 算出

• 當
  $n$ 非質數時：

根據 Bezout's Thm 有

$$ax+ny=1\Rightarrow a\cdot x\equiv 1\mathrm{mod}n$$
可用擴展歐幾里得演算法找到

$x=a^{-1}$

練習：

ZERO JUDGE a289 Modular Multiplicative Inverse

質數判斷

若數

$n>1$ 能被

$1$ 與

$n$ 以外的數整除，則

$n$ 為合數

bool is_p = true; // is_p := is it prime?
for (int i = 2; i <= sqrt(n); i++)
  if (n%i == 0) is_p = false;

練習：

CODEFORCES 1033B Square Difference

質因數分解

根據唯一分解定理，任何合數都能被分解成一些質數的積

算術基本定理

for (int p = 2; p <= sqrt(n); p++) {
  int t = 0;
  while (n%p == 0) n /= p, t++;
  if (t) factors.push_back({p, t});
}
if (n != 1) factors.push_back({n, 1});

練習：

CODEFORCES 1165D Almost All Divisors
CODEFORCES 1114C Trailing Loves (or L'oeufs?)
LeetCode 952 Largest Component Size by Common Factor

質數篩檢

數質數可以有效安撫緊張的情緒

Sieve of Eratosthenes

其精神是將質數的倍數都設為非質數

fill(is_p, is_p+maxn, true);
is_p[1] = false;

for (int n = 2; n < sqrt(maxn); n++) {
  if (!is_p[n]) continue;
  for (int m = n*n; m < maxn; m+=n) is_p[m] = false;
}

n*n 為對於所有

$x<n$，在數到

$n$ 之前

$x\cdot n$ 已經被篩過了。

記得檢查是否 n*n 溢位

練習：

UVa OJ 543 Goldbach’s Conjecture
UVa OJ 10140 Prime Distance
* ZERO JUDGE a007 判斷質數 ^[1]

線性篩法

fill(is_p, is_p+maxn, true);
is_p[1] = false;

for (int n = 2; n < maxn; n++) {
  if (is_p[n]) prime.push_back(n);
  
  for (int p: prime) {
    if (p*n >= maxn) break; // 超出篩檢範圍
    
    is_p[p*n] = false;
    if (n%p == 0) break;
  }
}

存在性

數到

$n$ 以前，它就已經被判定過是否為質數了。

根據唯一分解定理，

$n=p_1\cdot p_2\cdot ...\cdot p_k$，其中

$p_1\le p_2\le ...\le p_k$

令

$x_1=p_2\cdot ...\cdot p_k$，由於

$x_1<n$，所以必然有

$p_1$ 與

$x_1$ 相乘得到

$n$

$p_1$ 就是透過第二層迴圈得來的

就算

$p_1$ 整除

$x_1$， if (n%p == 0) break; 是之後才執行的，一定能湊得

$n$

唯一性

任意

$n$ 只會被判定過恰好一次

令

$x_i=\frac{n}{p_i}$，

$i\in \{1,2,..k\}$

$n$ 只會被數到

$x_1$ 時給判定，同樣用

$p_1$ 湊得

$n=p_1\cdot x_1$
對於

$x_i\ne x_1$， 在數到

$p_1$ 時

$x_i$ 就被

$p_1$ 給整除了，因而不會湊得

$p_i\cdot x_i$

$p_i$ 根據假設有

$p_i>p_1$

算術

對於大數字的運算，普通的做法不夠快，因此接下來將介紹快速的乘法及冪運算

Gauss's complex multiplication algorithm

對於複數

$a+bi$ 與

$c+di$ 相乘

$$(a+bi)\cdot (c+di)=(ac-bd)+(bc+ad)i$$
相較於

$4$ 次乘法

$ac,bd,bc,ad$

可只用

$3$ 次乘法

$\{\begin{array}{l}{k}_1=c(a+b)\\ k_2=b(c+b)\\ k_3=a(d-c)\end{array}$ 完成複數乘法：

$$\text{where }\{\begin{array}{l}\begin{array}{rlr}(ac-bd)& =ac+\mathbf{bc}-\mathbf{bc}-bd& =\underset{―}{c\cdot (a+b)}-\underset{―}{b\cdot (c+b)}=k_1-k_2\\ (bc+ad)& =\mathbf{ac}+bc+ad-\mathbf{ac}& =\underset{―}{c\cdot (a+b)}+\underset{―}{a\cdot (d-c)}=k_1+k_3\end{array}\end{array}$$
因此複數相乘

$$(a+bi)+(c+di)=(k_1-k_2)+(k_1+k_3)i$$

Karatsuba algorithm

對於整數

$\{\begin{array}{l}\begin{array}{rl}x& =am+b\\ y& =cm+d\end{array}\end{array}$ 相乘

受到複數乘法啟發，將整數

$x$ 分割成兩數

$a,b$

$$x\cdot y=ac\cdot m^2+(ad+bc)\cdot m+bd$$

能用

$3$ 次乘法

$\{\begin{array}{l}{z}_1=ac\\ z_2=bd\\ z_3=(a+b)(c+d)\end{array}$ 就完成運算：

$$\begin{array}{rl}(ad+bc)& =\mathbf{ac}+ad+bc+\mathbf{bd}-\mathbf{ac}-\mathbf{bd}\\ & =\underset{―}{(a+b)(c+d)}-\underset{―}{ac}-\underset{―}{bd}\\ & =z_3-z_1-z_2\end{array}$$

因此整數相乘

$$x\cdot y=z_1{m}^2+(z_3-z_1-z_2)m+z_2$$

範例

$\{\begin{array}{l}\begin{array}{rlrl}12345& =12& \cdot 1000& +345\\ 6789& =6& \cdot 1000& +789\end{array}\end{array}$相乘：

首先

$\{\begin{array}{ll}{z}_1=12\cdot 6& =72\\ z_2=345\cdot 789& =272205\\ z_3=(12+345)(6+789)=357\cdot 795& =283815\end{array}$

則

$\begin{array}{rl}12345\cdot 6789& =72\cdot {1000}^2+(283815-72-272205)\cdot 1000+272205\\ & =72000000+11538000+272205\\ & =83810205\end{array}$

似乎沒有比較快欸
那是因為你想的乘法，仍然是

$O(n^2)$ 的直式乘法

分治法

對於上述演算法，凡是遇到乘法運算，都使用同樣的演算法
並且對於數字的分割，總是分割成均等的兩半

例如上面

$6789$ 應分成

$67$ 和

$89$

$6789=67\cdot 100+89$

int k(int x, int y) {
  if(x < 10 || y < 10) return x*y;

  int len = min(log10(x), log10(y));
  int m = pow(10, len/2 + 1);

  auto [a, b] = div(x, m); // since c++17
  auto [c, d] = div(y, m);

  int z1 = k(a, c);
  int z2 = k(b, d);
  int z3 = k(a+b, c+d);

  return z1*m*m + (z3-z1-z2)*m + z2;
}

時間成本

$T(n)$ 約為

$$\begin{array}{rl}T(n)& =3\cdot T(\lceil n/2\rceil )+c_n\\ T(1)& =1+c_1\end{array}$$
此處

$n$ 為整數長度，

$c_k$ 為加法運算的成本
複雜度為

$O(3^{{\log }_{2}n})=O(n^{{\log }_{2}3})$

快速冪

求

$a$ 的

$n$ 次冪

$x=a^n$：

int x = 1;
while (n--) x *= a;

基於 D&C，若

$n$ 為偶數則

$a^n=a^{\frac{n}{2}}\times a^{\frac{n}{2}}$

• top-down

int fast_exp(int a, int n) {
  if (n == 1) return a;
  
  int x = fast_exp(a, n/2);
  return (x*x) * (n&1? a : 1); // 檢查是否為奇數
}

並且

$a^{n+m}=a^n\times a^m$

• bottom-up

int x = 1;
while (n) {
  if (n&1) x *= a;
  a *= a;
  n >>= 1;
}

複雜度從樸素的

$O(n)$ 改進至

$O({\log }_{2}n)$

練習：

UVa OJ 374 Big Mod
UVa OJ 10006 Carmichael Numbers
CODEFORCES 615D Multipliers

矩陣快速冪

對於方陣乘法，也能沿用快速冪的想法
例如求第

$n$ 項的費式數

$f_n$：

$${\left[\begin{array}{cc}1& 1\\ 1& 0\end{array}\right]}^n\cdot \left[\begin{array}{c}1\\ 0\end{array}\right]=\left[\begin{array}{c}{f}_{n+1}\\ f_n\end{array}\right]$$

int M[][2] = {{1, 1}, {1, 0}};
int f[] = {1, 0};

while (n) {
  if (n&1) {
    int t[] = {f[0], f[1]};
    f[0] = M[0][0] * t[0] + M[0][1] * t[1];
    f[1] = M[1][0] * t[0] + M[1][1] * t[1];
  }

  int p[][2] = {{M[0][0], M[0][1]}, {M[1][0], M[1][1]}};
  M[0][0] = p[0][0] * p[0][0] + p[0][1] * p[1][0];
  M[0][1] = p[0][0] * p[0][1] + p[0][1] * p[1][1];
  M[1][0] = p[1][0] * p[0][0] + p[1][1] * p[1][0];
  M[1][1] = p[1][0] * p[0][1] + p[1][1] * p[1][1];

  n >>= 1;
}

其複雜度

$O({\log }_{2}n)$

練習：

ZERO JUDGE b525 先別管這個了，你聽過turtlebee嗎？

---

1. 這根本不是基礎題= = ↩︎