結局一次不定方程式ax+by=cの整数解 | 高校数学の美しい物語を見れば導出できるのだけれど、具体的な導出結果が載っていなかったので自分用にメモ。
競技プログラミングでは正整数の場合しか出てこないと思うけれど、せっかくなので負の場合も考えている。
要約
0でない整数$a,,b$と整数$c$が与えられたとき、$ax + by = c$の整数解は以下のようになる。
$c$が$\mathrm{gcd}(a, b)$の倍数でないとき、$x,,y$は整数解を持たない。
$c$が$\mathrm{gcd}(a, b)$の倍数であるとき、$x,,y$は整数解を持ち、
\begin{align*}
x &= \frac{cx'_0 - bm}{\gcd(a, b)}\\
y &= \frac{cy'_0 + am}{\gcd(a, b)}\\
\end{align*}
ただし、$m$は任意の整数、$x'_0,,y'_0$は$ax' + by' = \gcd(a, b)$の特殊解。
$\gcd(a, b),,x'_0,,y'_0$は拡張互除法を用いることで一度に導出可能。
拡張互除法の計算量は$O(\log(\min(|a|, |b|)))$。
解説
$c$が$\mathrm{gcd}(a, b)$の倍数であるときのみ$x,,y$が整数解を持つ証明は一次不定方程式ax+by=cの整数解 | 高校数学の美しい物語を参照。
さて、拡張互除法は、0でない整数$a$と$b$が与えられたとき、次の値を一度に求めるアルゴリズムである。
- $a$と$b$を共に割り切れる整数のうち最大の正整数（これを最大公約数と呼び通常$\gcd(a, b)$と表す）
- $ax' + by' = \gcd(a, b)$を満たす整数解$x',,y'$（無数にある）のうちのどれか一つの組（これを特殊解と呼ぶ）
$ax' + by' = \gcd(a, b)$の特殊解を$x'_0,,y'_0$と呼ぶことにしよう。
拡張互除法の具体的なアルゴリズムは後で説明することにして、ここでは$\gcd(a, b),, x'_0,,y'_0$の具体的な値が求まったものとする。
$c$が$\gcd(a, b)$の倍数であるとき、$ax + by = c$の特殊解$x_0,,y_0$も$x'_0,,y'_0$を使って表すことができる。
すなわち、
\begin{align*}
ax'_0 + by'_0 &= \gcd(a, b)\\
\frac{c}{\gcd(a, b)}(ax'_0 + by'_0) &= \frac{c}{\gcd(a, b)}\gcd(a, b)\\
a \frac{c}{\gcd(a, b)} x'_0 + b \frac{c}{\gcd(a, b)} y'_0 &= c\\
x_0 &= \frac{c}{\gcd(a, b)}x'_0\\
y_0 &= \frac{c}{\gcd(a, b)}y'_0\\
\end{align*}
特殊解$x_0,,y_0$を用いると、一般解$x,,y$を求めることができる。
すなわち、
\begin{align*}
ax + by &= c\\
ax_0 + by_0 &= c\\
(ax + by) - (ax_0 + by_0) &= c - c\\
a(x - x_0) + b(y - y_0) &= 0\\
a(x - x_0) &= -b(y - y_0)\\
\frac{a}{\gcd(a, b)}(x - x_0) &= -\frac{b}{\gcd(a, b)}(y - y_0)\\
\end{align*}
$\frac{a}{\gcd(a, b)}$と$-\frac{b}{\gcd(a, b)}$は互いに素なので、両辺は任意の整数$m$を用いて次のような形になるときしか等しくならない。
\begin{align*}
\frac{a}{\gcd(a, b)}\left(-\frac{b}{\gcd(a, b)}m\right) &= -\frac{b}{\gcd(a, b)}\left(\frac{a}{\gcd(a, b)}m\right)\\
\end{align*}
すなわち、
\begin{align*}
x - x_0 &= -\frac{b}{\gcd(a, b)}m\\
y - y_0 &= \frac{a}{\gcd(a, b)}m\\
x &= x_0 - \frac{b}{\gcd(a, b)}m\\
y &= y_0 + \frac{a}{\gcd(a, b)}m\\
x &= \frac{c}{\gcd(a, b)}x'_0 - \frac{b}{\gcd(a, b)}m\\
y &= \frac{c}{\gcd(a, b)}y'_0 + \frac{a}{\gcd(a, b)}m\\
x &= \frac{cx'_0 - bm}{\gcd(a, b)}\\
y &= \frac{cy'_0 + am}{\gcd(a, b)}\\
\end{align*}
拡張互除法
以下のコードはC++11以降で正しく動く。
#include <iostream>
#include <tuple>
// 返り値の0番目の要素が最大公約数、1番目と2番目の要素が特殊解
std::tuple<long long, long long, long long> extgcd(long long prev_r, long long r) {
long long prev_s = 1;
long long prev_t = 0;
long long s = 0;
long long t = 1;
while (r != 0) {
const long long q = prev_r > 0 ? prev_r / r : r > 0 ? -((-prev_r - 1 + r) / r) : (-prev_r - 1 + -r) / -r;
const long long next_r = prev_r - q * r;
prev_r = r;
r = next_r;
const long long next_s = prev_s - q * s;
prev_s = s;
s = next_s;
const long long next_t = prev_t - q * t;
prev_t = t;
t = next_t;
}
return std::make_tuple(prev_r, prev_s, prev_t);
}
int main() {
const long long a = -33447ll;
const long long b = 90629ll;
auto result = extgcd(a, b);
std::cout << a << " * " << std::get<1>(result) << " + " << b << " * " << std::get<2>(result) << " = " << std::get<0>(result) << std::endl;
// -33447 * -7400 + 90629 * -2731 = 1
return 0;
}
解説
まず通常の互除法の復習から行こう。
通常の互除法のうち我々がよく知るものは、正整数$a,,b$が与えられたとき$\gcd(a, b)$のみしか求めないもので、以下のようにして計算する。（$x,\texttt{%},y$は$x$を$y$で割った余りとする）
\begin{alignat*}{3}
r_0 &= a &&\\
r_1 &= b &&\\
r_2 &= r_0 \,\texttt{%}\, r_1 &&\\
r_3 &= r_1 \,\texttt{%}\, r_2 &&\\
&\ \ \vdots &&\\
r_{i + 1} &= r_{i - 1} \,\texttt{%}\, r_i &&\\
&\ \ \vdots &&\\
r_k &= r_{k - 2} \,\texttt{%}\, r_{k - 1} &&= \gcd(a, b)\\
r_{k + 1} &= r_{k - 1} \,\texttt{%}\, r_k &&= 0
\end{alignat*}
$a < b$であっても$a ,\texttt{%}, b = a$なので、$a$と$b$が入れ替わってうまく動いてくれることに注意。
ここから、0でない整数$a,,b$が与えられたとき$\gcd(a, b)$を求め、さらに$ax' + by' = \gcd(a, b)$の特殊解$x'_0,,y'_0$も求めるアルゴリズムを考える。
このままだと扱いづらいので式を変形する。
\begin{alignat*}{3}
r_0 &= q_1 r_1 + r_2 && (\text{ただし}\,0 \leq r_2 < r_1)\\
r_1 &= q_2 r_2 + r_3 && (\text{ただし}\,0 \leq r_3 < r_2)\\
&\ \ \vdots &&\\
r_{k - 1} &= q_k r_k + r_{k + 1} && (\text{ただし}\,0 \leq r_{k + 1} < r_k)\\
\end{alignat*}
$q_i$は商を表す整数であり、よくある商と余りの関係を表す式になっている。
この辺りで$a$と$b$が0でない整数の場合にも式を適用できるように拡張しておきたい。
一般に、$a$と$b$が0でない整数のとき、$a=qb+r,(\text{ただし},0 \leq r < b)$という条件では整数$q,,r$を決めることはできない。
そこで、ここでは最小非負剰余という定義を用いて$q$と$r$が一意に定まるようにする。
最小非負剰余とは$a=qb+r,(\text{ただし},0 \leq r < |b|)$で決まる$r$のことである。
$r$が決まると同時に$q$も定まり、$a>0$では$\frac{a}{b}$を0方向へ丸めた値、$a<0$では$\frac{a}{b}$を無限大方向へ丸めた値となる。
すなわち、
\begin{alignat*}{3}
r_0 &= q_1 r_1 + r_2 && (\text{ただし}\,0 \leq r_2 < |r_1|)\\
r_1 &= q_2 r_2 + r_3 && (\text{ただし}\,0 \leq r_3 < |r_2|)\\
&\ \ \vdots &&\\
r_{k - 1} &= q_k r_k + r_{k + 1} && (\text{ただし}\,0 \leq r_{k + 1} < |r_k|)\\
\end{alignat*}
さて、これらの関係は行列の積として表すことができる。
すなわち、
\begin{align*}
\left(\begin{array}{c}
r_0\\
r_1
\end{array}\right) &= \left(\begin{array}{cc}
q_1 & 1\\
1 & 0
\end{array}\right)\left(\begin{array}{c}
r_1\\
r_2
\end{array}\right)\\
\left(\begin{array}{c}
r_1\\
r_2
\end{array}\right) &= \left(\begin{array}{cc}
q_2 & 1\\
1 & 0
\end{array}\right)\left(\begin{array}{c}
r_2\\
r_3
\end{array}\right)\\
&\ \ \vdots\\
\left(\begin{array}{c}
r_{k - 1}\\
r_k
\end{array}\right) &= \left(\begin{array}{cc}
q_k & 1\\
1 & 0
\end{array}\right)\left(\begin{array}{c}
r_k\\
r_{k + 1}
\end{array}\right)\\
\end{align*}
ということは、
\begin{align*}
\left(\begin{array}{c}
r_0\\
r_1
\end{array}\right) &= \left(\begin{array}{cc}
q_1 & 1\\
1 & 0
\end{array}\right)\left(\begin{array}{cc}
q_2 & 1\\
1 & 0
\end{array}\right)\cdots\left(\begin{array}{cc}
q_k & 1\\
1 & 0
\end{array}\right)\left(\begin{array}{c}
r_k\\
r_{k + 1}
\end{array}\right)\\
\left(\begin{array}{c}
a\\
b
\end{array}\right) &= \left(\begin{array}{cc}
q_1 & 1\\
1 & 0
\end{array}\right)\left(\begin{array}{cc}
q_2 & 1\\
1 & 0
\end{array}\right)\cdots\left(\begin{array}{cc}
q_k & 1\\
1 & 0
\end{array}\right)\left(\begin{array}{c}
\gcd(a, b)\\
0
\end{array}\right)\\
\end{align*}
ここで、
\left|\begin{array}{cc}
q_i & 1\\
1 & 0
\end{array}\right| \neq 0
なので逆行列が存在し、
\begin{align*}
\left(\begin{array}{cc}
q_k & 1\\
1 & 0
\end{array}\right)^{-1}\cdots\left(\begin{array}{cc}
q_2 & 1\\
1 & 0
\end{array}\right)^{-1}\left(\begin{array}{cc}
q_1 & 1\\
1 & 0
\end{array}\right)^{-1}\left(\begin{array}{c}
a\\
b
\end{array}\right) &= \left(\begin{array}{c}
\gcd(a, b)\\
0
\end{array}\right)\\
\left(\begin{array}{cc}
0 & 1\\
1 & -q_k
\end{array}\right)\cdots\left(\begin{array}{cc}
0 & 1\\
1 & -q_2
\end{array}\right)\left(\begin{array}{cc}
0 & 1\\
1 & -q_1
\end{array}\right)\left(\begin{array}{c}
a\\
b
\end{array}\right) &= \left(\begin{array}{c}
\gcd(a, b)\\
0
\end{array}\right)\\
\end{align*}
ということが言える。
さて、ここで、
\begin{align*}
A_1 &= \left(\begin{array}{cc}
1 & 0\\
0 & 1
\end{array}\right)\\
A_{i + 1} &= \left(\begin{array}{cc}
0 & 1\\
1 & -q_i
\end{array}\right) A_i
\end{align*}
として行列$A_i$を定義する。
つまり、
\begin{align*}
A_{k+1}\left(\begin{array}{c}
a\\
b
\end{array}\right) &= \left(\begin{array}{c}
\gcd(a, b)\\
0
\end{array}\right)\\
\end{align*}
である。
また、
\begin{align*}
\left(\begin{array}{cc}
u_i & v_i\\
s_i & t_i
\end{array}\right) &= A_i
\end{align*}
として整数$u_i,,v_i,,s_i,,t_i$を定義する。
すると、
\begin{align*}
A_{i + 1} &= \left(\begin{array}{cc}
0 & 1\\
1 & -q_i
\end{array}\right)\left(\begin{array}{cc}
u_i & v_i\\
s_i & t_i
\end{array}\right)\\
&= \left(\begin{array}{cc}
s_i & t_i\\
u_i - q_i s_i & v_i - q_i t_i
\end{array}\right)
\end{align*}
なので、
\begin{align*}
u_{i + 1} &= s_i\\
v_{i + 1} &= t_i\\
s_{i + 1} &= u_i - q_i s_i\\
t_{i + 1} &= v_i - q_i t_i
\end{align*}
より、
\begin{align*}
u_i &= s_{i - 1}\\
v_i &= t_{i - 1}\\
s_{i + 1} &= s_{i - 1} - q_i s_i\\
t_{i + 1} &= t_{i - 1} - q_i t_i
\end{align*}
が言える。
ということは、
\begin{align*}
A_i &= \left(\begin{array}{cc}
s_{i - 1} & t_{i - 1}\\
s_i & t_i
\end{array}\right)
\end{align*}
であるから、
\begin{align*}
\left(\begin{array}{cc}
s_k & t_k\\
s_{k + 1} & t_{k + 1}
\end{array}\right)\left(\begin{array}{c}
a\\
b
\end{array}\right) &= \left(\begin{array}{c}
\gcd(a, b)\\
0
\end{array}\right)\\
\end{align*}
より、
\begin{align*}
a s_k + b t_k &= \gcd(a, b)
\end{align*}
であり、$s_k,,t_k$が$ax' + by' = \gcd(a, b)$の特殊解になっていることがわかる。
つまり、これまでの考察をまとめると、
\begin{align*}
r_0 &= a\\
s_0 &= 1\\
t_0 &= 0\\
r_1 &= b\\
s_1 &= 0\\
t_1 &= 1\\
q_1 &= \left\{\begin{array}{ll}
\cfrac{r_0}{r_1}\text{を}0\text{方向に丸めた値} & \text{if}\ r_0 > 0\\
\cfrac{r_0}{r_1}\text{を無限大方向に丸めた値} & \text{if}\ r_0 < 0\\
\end{array}\right.\\
r_2 &= r_0 - q_1 r_1\\
s_2 &= s_0 - q_1 s_1\\
t_2 &= t_0 - q_1 t_1\\
&\ \ \vdots\\
q_{i - 1} &= \left\{\begin{array}{ll}
\cfrac{r_{i - 2}}{r_{i - 1}}\text{を}0\text{方向に丸めた値} & \text{if}\ r_{i - 2} > 0\\
\cfrac{r_{i - 2}}{r_{i - 1}}\text{を無限大方向に丸めた値} & \text{if}\ r_{i - 2} < 0\\
\end{array}\right.\\
r_i &= r_{i - 2} - q_{i - 1} r_{i - 1}\\
s_i &= s_{i - 2} - q_{i - 1} s_{i - 1}\\
t_i &= t_{i - 2} - q_{i - 1} t_{i - 1}\\
&\ \ \vdots\\
q_{k - 1} &= \left\{\begin{array}{ll}
\cfrac{r_{k - 2}}{r_{k - 1}}\text{を}0\text{方向に丸めた値} & \text{if}\ r_{k - 2} > 0\\
\cfrac{r_{k - 2}}{r_{k - 1}}\text{を無限大方向に丸めた値} & \text{if}\ r_{k - 2} < 0\\
\end{array}\right.\\
r_k &= r_{k - 2} - q_{k - 1} r_{k - 1} = \gcd(a, b)\\
s_k &= s_{k - 2} - q_{k - 1} s_{k - 1} = x'_0\\
t_k &= t_{k - 2} - q_{k - 1} t_{k - 1} = y'_0\\
q_k &= \left\{\begin{array}{ll}
\cfrac{r_{k - 1}}{r_k}\text{を}0\text{方向に丸めた値} & \text{if}\ r_{k - 1} > 0\\
\cfrac{r_{k - 1}}{r_k}\text{を無限大方向に丸めた値} & \text{if}\ r_{k - 1} < 0\\
\end{array}\right.\\
r_{k + 1} &= r_{k - 1} - q_k r_k = 0\\
\end{align*}
である。
これをコードに落とすと上で書いたようなものになる。
感想
競技プログラミング、きちんと証明しようと思うととんでもない量になる定理が出がちでは？
競プロ上位者はだいたい数オリ上位者も兼ねているという話も聞きますし、この程度で音を上げるようではダメなんでしょうね。