• 本文總計約 8000 字，閱讀大約需要 30 分鐘。
• 警告！警告！警告！本文有大量的 \(\LaTeX\) 公式渲染，可能會導致載入異常緩慢！

前言

擴充套件 Euclid 演算法是數論中最經典，也是最簡單（？）的演算法了，但是它的應用又十分的廣泛。所以介紹它是十分有必要的。

數論的知識真的是非常奇怪呢（bushi），因為它們感性地理解起來非常輕鬆，但是要嚴謹證明就很難 QAQ，筆者在文中加入了大量的數學證明，儘管可能有讀者不愛看因為其實我也不愛看，但是基於數學嚴謹的思想，我還是把它們加了進來，這可能會導致瀏覽器渲染起來蠻麻煩的，不過還是請讀者多多包涵了。

題目引入

給定三個整數 \(a,b,c\)

，求關於 \(x,y\) 的二元一次不定方程 \(ax+by=c\) 的一組整數解，注意該方程可能沒有整數解。

樣例輸入 \(1\)：\(a=3,b=2,c=7\)；
樣例輸出 \(1\)：\(x=1,y=2\)。注意可能有其它的滿足要求的解，如 \(x=3,y=-1\)。

樣例輸入 \(2\)：\(a=2,b=-4,c=1\)；
樣例輸出 \(2\)：無解，因為 \(2x-4y\) 永遠是偶數，但 \(c=1\)，是奇數。

最大公約數及 Euclid 演算法

最大公約數的概念及性質

即使最大公約數在我們小學的時候就已經接觸了，但我們還是需要重新介紹一下它：

如果對於兩個整數 \(a,b\)，若 \(m\)

是最大的既可以整除 \(a\)（記作 \(m\mid a\)），又可以整除 \(b\) 的正整數，即 \(m\) 是 \(a,b\) 最大的公約數，那麼我們稱 \(m\) 為 \(a,b\) 的最大公約數（好像很廢話的概念），記作 \(m=\gcd(a,b)\)。

特別地，若兩個整數 \(a,b\) 滿足 \(\gcd(a,b)=1\)，則我們稱 \(a\) 與 \(b\) 互質，記作 \(a\perp b\)。

顯然 \(\gcd\) 函式有如下的性質，此處不加以證明：

• 交換律：\(\gcd(a,b)=\gcd(b,a)\)，特別地，\(\gcd(a,0)=\gcd(0,a)=1\)；
• 可以提出常數：\(\gcd(ka,kb)=k\gcd(a,b)\ (k\in\mathbb{Z}\)
  )；
• 對加減法及取模自由：\(\gcd(a,b)=\gcd(a,b\pm a)=\gcd(a,b\pm ka)=\gcd(a,b \bmod a)\ (k\in\mathbb{Z})\)。

Euclid 演算法的內容

現在我們已經知道了最大公約數的基本概念了，那麼我們怎麼求任意兩個給定的整數 \(m,n(m\ge n)\) 的最大公約數呢？

答曰：我會 \(\mathcal{O}(n)\) 暴力列舉！

然而，假如我們令 \(n\le 10^9\) 甚至是 \(n\le 10^{18}\)，那麼 \(\mathcal{O}(n)\) 的列舉是一定會 TLE 的 QwQ，我們就需要一種更快的演算法來求出上述問題中的 \(\gcd(m,n)\)。於是，Euclid 演算法應運而生。

Euclid 演算法，又稱輾轉相除法，顧名思義，是由古希臘數學家 \(\mathcal{Euclid}\) 發明的演算法，其步驟如下：

• 對於兩個整數 \(m,n\)，如果有 \(n\ne0\)，那麼我們求出 \(r=m\bmod n\) 的值；
• 接下來，使 \(m=n,n=r\)；
• 重複上述兩個步驟，直到 \(n=0\) 時，\(m\) 的值就是我們所求的最大公約數。

例如，我們要求 \(\gcd(1081,897)\)，步驟如下：

• \(n=897\ne0\)，故我們令 \(m=897\)，\(n=1081\bmod897=184\)；原問題即轉化為求 \(\gcd(897,184)\)；
• \(n=184\ne0\)，故我們令 \(m=184\)，\(n=897\bmod184=161\)；原問題即轉化為求 \(\gcd(184,161)\)；
• \(n=161\ne0\)，故我們令 \(m=161\)，\(n=184\bmod161=23\)；原問題即轉化為求 \(\gcd(161,23)\)；
• \(n=23\ne0\)，故我們令 \(m=23\)，\(n=161\bmod23=0\)；此時 \(n=0\)，所以 \(\gcd(1081,897)=m=23\)。

程式碼

Euclid 演算法的程式碼非常好寫，這裡給出迭代求最大公約數的程式碼：

#include <cstdio>
#include <algorithm>  //有用於求 gcd 的標準函式 __gcd(int, int)
using namespace std;

inline int Euclid(int x, int y) {  //輾轉相除法求 gcd(x, y)，inline 防止無效資訊入棧
	int r;
	
	while(y) {
		r = x % y;
		x = y;
		y = r;
	}
	
	return x;
}

long long m, n;

int main(void) {
	scanf("%lld%lld", &m, &n);
	printf("%lld", Euclid(m, n));
	return 0;
}
//by CaO

當然，C++ 中 <algorithm> 庫裡有 __gcd(int, int) 函式，功能與上述函式相同，但我不知道 NOIP 中讓不讓用 QwQ。

時間複雜度分析

數學證明警告！
有一種非常經典的證明其時間複雜度的方法，如下：

我們構造數列 \(\{U\}\)，其通項公式為：

\[U_i=\begin{cases} m & i=0\\ n & i=1\\ U_{i-2}\bmod U{i-1} & \text{otherwise.} \end{cases} \]

定義 Fibonacci 數列 \(\{F\}\) 通項公式為：

\[F_i=\begin{cases} 1 & i\le 1 \\ F_{i-1}+F_{i-2} & i>2 \end{cases}\ (i=0,1,\cdots) \]

顯然數列 \(\{U\}\) 中存在一項 \(U_x=0\)。

故有 \(F_0=U_x=0,F_1\ge U_x\)。

又因為 \(U_i \bmod U_{i+1}=U_{i+2}\)，所以 \(U_{i}\ge U_{i+1}+U_{i+2}\)。

利用數學歸納法可得：\(U_i\ge F_{x-i+1}\)，故 \(n=U_1\ge F_{x}\)。

注意到 \(F_x>\dfrac{(\phi+1)^x}{\sqrt{5}}\)，其中 \(\phi=\dfrac{\sqrt{5}-1}{2}\)，所以 \(x<\log_{(\phi+1)}(\sqrt{5}n)\)。

故 Euclid 演算法迭代的次數不超過 \(\log_{(\phi+1)}(\sqrt{5}n)\) 次，時間複雜度為 \(\mathcal{O}(\log n)\)。\(\Box\)

這就是說，Euclid 演算法相比較於 \(\mathcal{O}(n)\) 的暴力列舉，有了很大程度的優化。

不定方程概念的引入

不定方程

一個含有多個未知數的整係數方程稱為不定方程，又稱 Diophantine 方程。即，形如：\(\sum_{i=1}^{n}a_ix_i^{p_i}=b\) 的方程即為丟番圖方程，其中 \(a_i,p_i,b\) 均為整數（\(i=1,2,\cdots,n\)）。如果這個方程只有兩個未知數，且未知數的次數均為 \(1\)，那麼我們稱其為二元一次不定方程。

注意，本文在討論不定方程時，預設指的是二元一次不定方程，即形如 \(ax+by=c\) 的方程。

不定方程有解問題 — — Bézout 定理

正如題目引入中所說的，二元一次不定方程是不一定有解的，例如 \(2x-4y=1\) 就沒有整數解。所以，我們應該如何判斷，一個不定方程是否有解呢？

針對這個問題，早在 \(17\) 世紀，法國數學家 \(\mathcal{B\acute{e}zout}\) 就給出了結論：

關於 \(x,y\) 的不定方程 \(ax+by=c\) 有整數解，當且僅當 \(\gcd(a,b)\mid c\)。或者用另一種表述：\(ax+by=\gcd(a,b)\) 一定有整數解。這就是著名的 Bézout 定理，或者國內又譯為裴蜀定理或貝祖定理。證明如下：（數學證明警告！）

先證充分性，即若 \(\gcd(a,b) \mid c\)，則 \(ax+by=c\) 有整數解。

設 \(S\) 是最小的能用 \(ax+by\) 表示的正整數，則有
\(a\bmod S=a-S\left \lfloor \dfrac{a}{S} \right \rfloor=a-(ax+by)\left \lfloor \dfrac{a}{S} \right \rfloor=a \left(1-\left \lfloor \dfrac{a}{S} \right \rfloor \right)-b\left(\left \lfloor \dfrac{a}{S} \right \rfloor\right)\le S\)。

又因為 \(S\) 是最小的能用 \(ax+by\) 表示的正整數，所以 \(a\bmod S=0\)，即 \(S\mid a\)；

同理 \(S\mid b\)。即 \(S\) 是 \(a,b\) 的公約數，由最大公約數定義得 \(S\le \gcd(a,b)\)。

又因為 \(\gcd(a,b)\mid (ax+by)\)，所以 \(\gcd(a,b)\mid S\)，故 \(\gcd(a,b) \le S\)。

綜上 \(S=\gcd(a,b)\)，即 \(ax+by=\gcd(a,b)\) 有整數解。

且對於任意 \(c=kS\ (k\in\mathbb{Z})\)，都存在 \(x'=kx,y'=ky\)，使得 \(ax'+by'=c\)，即 \(ax+by=c\) 存在正整數解。

再證必要性，即若 \(ax+by=c\) 有整數解，則 \(\gcd(a,b) \mid c\)。

假設 \(x_0,y_0\) 是 \(ax+by=\gcd(a,b)\) 的整數解，則任意的 \(ax+by=c\) 的整數解 \(x,y\)，都有 \(x=\dfrac{c}{\gcd(a,b)}x_0,y=\dfrac{c}{\gcd(a,b)}y_0\)。

若 \(\gcd(a,b)\nmid c\)，則 \(x,y\notin\mathbb{Z}\)，與 \(x,y\) 為整數矛盾。

綜上所述，不定方程 \(ax+by=c\) 有整數解當且僅當 \(\gcd(a,b)\mid c\)。\(\Box\)

Bézout 定理的推廣

當然，即使 Bézout 定理是關於二元一次不定方程的，但我們依舊可以將其推廣到 \(n\) 元的不定方程：

關於 \(x_1,x_2,\cdots,x_n\) 的 \(n\) 元一次不定方程 \(\sum_{i=1}^{n}a_ix_i=b\) 有整數解，當且僅當 \(\gcd_{i=1}^{n}(a_i)\mid b\)。證明留給讀者作為練習。

不定方程的解

擴充套件 Euclid 演算法的引入

對於不定方程 \(ax+by=c\)，我們可以通過 Euclid 演算法求出 \(\gcd(a,b)\) 的值，然後用 Bézout 定理判斷其是否存在整數解。

但是，我們只是知道了一個不定方程的可解性，又如何求它的具體解呢？

其實，這個問題通過對 Euclid 演算法進行一些小小的拓展、變形，就可以求出來。也就是接下來的擴充套件 Euclid 演算法。

擴充套件 Euclid 演算法的內容

我們現在考慮最簡單的一種情況：\(ax+by=\gcd(a,b)\)，應該如何求其整數解呢？

首先，小學的數學知識就告訴我們，商乘除數加餘數等於被除數，這一點是永遠成立的，所以我們用 \(b\) 除以 \(a\)，有 \(b=a\left \lfloor \dfrac{b}{a} \right \rfloor + (b\bmod a)\)。代入原方程，有 \(ax+\left(a\left \lfloor \dfrac{b}{a} \right \rfloor + (b\bmod a)\right)y=\gcd(a,b)\)。

重新合併，得到了 \(a\left(\left \lfloor \dfrac{b}{a} \right \rfloor y+x\right)+(b\bmod a)y=\gcd(a,b)\)。

右邊由 \(\gcd\) 的性質有 \(\gcd(a,b)=\gcd(a,a\bmod b)\)，我們令左邊 \(x'=y,y'=\left \lfloor \dfrac{b}{a} \right \rfloor y+x\)，則有 \((a\bmod b)x'+ay'=\gcd(a,a\bmod b)\)，此時，我們就得到了一個新的，係數範圍更小的不定方程。

如果我們求出了方程中的 \(x',y'\)，就可以代入回去，得到 \(y=x',x=y'-\left \lfloor \dfrac{b}{a} \right \rfloor x'\)。整個演算法就可以遞迴地進行了。

遞迴的出口為 \(b=0\)，此時我們有 \(x=1,y=0\)。

像這樣，我們利用了係數與係數之間 “輾轉相除” 的方式，求出了 \(ax+by=\gcd(a,b)\) 的一組解，同時，我們也可以求出 \(\gcd(a,b)\) 的值，即遞迴出口處 \(a\) 的值。

這種演算法是通過 Euclid 演算法變形而來，所以我們稱其為擴充套件 Euclid 演算法，簡稱擴歐。

這種演算法只能求出來不定方程的一組解，那麼我們怎麼求出其它的解呢？

假如我們通過擴歐求出的解為 \(x_0,y_0\)，則 \(x=x_0+k\dfrac{bc}{\gcd(a,b)},y=y_0-k\dfrac{ac}{\gcd(a,b)}\ (k\in \mathbb{Z})\) 也是不定方程的一組解。同時，我們可以找到其中 \(x\) 的最小正值為 \(\left(x\bmod \dfrac{bc}{\gcd(a,b)}+\dfrac{bc}{\gcd(a,b)}\right)\bmod \dfrac{bc}{\gcd(a,b)}\)。\(y\) 值同理。

程式碼

#include <cstdio>
using namespace std;

int a, b, c, k, x, y;
int exgcd(int a, int b, int& x, int& y) {  //擴歐求不定方程
	if(b == 0) {
		x = 1;
		y = 0;
		return a;  //遞迴出口
	}
	
	int ans = (b, a % b, y, x);
	y -= (a / b * x);  //狀態轉移
	return ans;
}

int main(void) {
	scanf("%d%d%d", &a, &b, &c);
	
	if(c % exgcd(a, b, x, y) != 0) {  //根據裴蜀定理，此時方程無整數解
		printf("No Solution.\n");
	}
	
	else {
		k = c / exgcd(a, b, x, y);
		printf("x=%d, y=%d\n", k * x, k * y);
	}
	
	return 0;
}
//by CaO

時間複雜度分析

擴歐的演算法與 Euclid 演算法步驟基本相同，時間複雜度依舊為 \(\mathcal{O}(\log n)\)。

擴充套件 Euclid 演算法的應用

擴歐演算法主要用在同餘方程之中：我們可以將一個類似 \(x\equiv a\pmod{n}\) 的方程寫為 \(x-nk=a\) 的形式，然後對 \(x,k\) 用擴歐解出一組解，找到滿足條件的 \(x\)。

也就可以應用於求模運算意義下的乘法逆元和中國剩餘定理（CRT 問題）的應用了。

例題

• 洛谷 P5656【模板】二元一次不定方程 (exgcd)
  這是擴歐的板子題，只是要特殊判斷一下輸出的型別還挺麻煩的 QAQ。
• 洛谷 P1516 青蛙的約會
  在模運算意義下的擴歐解不定方程，其實也很水（雖然是道藍題）。