edo1z blog

プログラミングなどに関するブログです

C++ 拡張ユークリッドの互除法

ユークリッドの互除法は2つの自然数の最大公約数を求める際に利用します。

  • a % b == r の場合、a, bの最大公約数は、b, rの最大公約数と等しいです。
  • b % r == r2、r % r2 == r3と繰り返すことで、最終的に余りが0になったときの割ってる方の数が最大公約数になります。

サンプルコード

#include <bits/stdc++.h>
using namespace std;
using LL = long long;

LL gcd(LL a, LL b) {
  if (b == 0) return a;
  return gcd(b, a % b);
}

int main() {
  cout << gcd(1029, 1071) << endl;
}

結果

21

拡張ユークリッドの互除法

参考

1次不定方程式を解く

  • 一次不定方程式 ax + by = gcd(a, b)の(x, y)を出すときに、ユークリッドの互除法が使えます。
  • 一次不定方程式についての説明は、ここに説明があります。

解き方

  • a % b = rの場合、gcd(a, b) == gcd(b, r)です。
  • ということは、ax + by = gcd(a,b) = bx' + ry' です。
  • ということは、下記のようになっていきます。ユークリッド互除法と同様です。
ax     + by     = gcd(a, b) ... (1)
bx'    + ry'    = gcd(b, r) ... (2)
rx''   + r'y''  = gcd(r, r')
r'x''' + 0y'''  = gcd(r', 0)
  • 最後のx'''は出せます。
  • r'x''' + 0y''' = gcd(r', 0) = r'
  • つまり、x''' == 1です。
  • y'''は、何でもいいです。
  • x'とy'から1つ上のxとyが出せれば、再帰的に最後のx'''とy'''から、最初のxとyが出せます。
  • rは、a % bで、a % bは、a - (a / b) * b です。
    • c++だと、整数同士の割り算は少数点以下を切り捨てるので、上記式が成り立ちます。
  • 上記を、(2)の式に代入すると、 下記のようになります。
bx'    + ry'    = gcd(b, r) ... (2)
= bx' + (a -(a/b)b)y'
= bx' + ay' - (a/b)by'
= ay' + b(x' - (a/b)y')
  • 上記と、1つ上の式(1)を比べると、下記になります。
x == y'
y == x' - (a/b)y'
  • これで再帰的に解くことができます。

コードサンプル

#include <bits/stdc++.h>
using namespace std;
using LL = long long;

LL extgcd(LL a, LL b, LL &x, LL &y) {
  if (b == 0) {
    x = 1;
    y = 0;
    return a;
  }
  LL gcd = extgcd(b, a % b, x, y);
  LL oldX = x;
  x = y;
  y = oldX - a / b * y;
  return gcd;
}

int main() {
  LL a = 1029, b = 1071, x, y;
  LL gcd = extgcd(a, b, x, y);
  printf("%d x %d + %d x %d = %d", a, x, b, y, gcd);
}

結果

1029 x 25 + 1071 x -24 = 21