$x,y$ の最大公約数 $\mathrm{gcd}(x,y)$ を計算します（最大公約数は英語でgreatest common divisorですので、gcdと略されます）。
なお、$x=y=0$ の場合には、任意の整数が $x,y$ の公約数となります。そのため公約数の最大値は存在せず、$\mathrm{gcd}(0,0)$は通常定義しません（イデアル論的な観点から、$0$と定義することもあります。）
数学的背景
$d$ が $n$ を割り切るとき、$d\mid n$ と書きます。次が成り立ちます：
\[
d\mid x \text{かつ} d\mid y \iff d\mid x \text{かつ} d\mid x – ky
\]
これは $d$ の倍数同士の和・差が $d$ の倍数であることから、簡単に証明できます。さらにこの定理から、$\mathrm{gcd}$ の性質\[\mathrm{gcd}(x,y) = \mathrm{gcd}(x, y-kx)\]が成り立つことが分かります。このことを用いると、$2$ 数の最大公約数を次のように計算できます：
\[
\mathrm{gcd}(345,100) = \mathrm{gcd}(100,45) = \mathrm{gcd}(45,10) = \mathrm{gcd}(10,5) = \mathrm{gcd}(5,0) = 5.
\] この計算方法を、Euclidの互除法といいます。大きい数を、小さい数で割った剰余に置き換えることを繰り返しています。
計算量
$a>b>0$ であり，$a$ を $b$ で割った余りが $c$ であるとします。
\[
a \geq 2b \implies b+c \leq b+b \leq \frac23(a+b),\qquad
a \leq 2b \implies b+c \leq b+(a-b) \leq \frac23(a+b)
\]より、Euclidの互除法の計算を1ステップ進めると、2数の和が $\frac23$ 倍以下になることが分かります。したがって、$\gcd(x,y)$ の計算は $O(\log(x+y))$ で行えます。1回目の除算でどちらも $\min(x,y)$ 以下になるため、$O(\log \min(x,y))$ とすることも可能です。（上の証明から、計算回数の上界が $\log_{1.5}(x+y)$ と評価できますが，この底 $1.5$ は精密化すると，黄金比$\frac{1+\sqrt{5}}{2} = 1.618$ まで改善できます。）
我々の見慣れている $2, 3$ 桁の整数であれば、素因数から最大公約数を求めるのが簡潔そうに見えますが、桁数が大きくなるとEuclidの互除法の方が圧倒的に性能が高く、計算機では数百桁といった数のgcdですら瞬時に計算することができます。
実装
gcd関数は、math・numpyなどのライブラリにも備わっていますが、ここでは自作してみます。
def gcd(a, b):
while b:
a, b = b, a % b
return adef gcd(a, b):
return gcd(b, a % b) if b else aループによる実装、再帰による実装の2つを挙げていますが、どちらも処理内容は同じです。
なお、関数呼び出しのオーバーヘッドの分、ループによる実装の方がパフォーマンスは良いです。
（関数の1回あたりの処理がある程度重いならば、呼び出しの負荷は相対的に無視量になりますが、今回は1回除算をする程度ですので、関数の呼び出しの差が相対的に多くなります。パフォーマンスについては別記事で、計測結果などを扱っていく予定です。）
関連：最小公倍数
$2$ 数 $a,b$ の最小公倍数を、$\mathrm{lcm}(a,b)$ と書きます（least common multiple）。$\mathrm{gcd}$ との間に次の関係があります：
\[ \mathrm{lcm}(a,b) = \frac{ab}{\mathrm{gcd}(a,b)}.
\] このことは、素因数分解を考えることで証明されます（3数以上の最大公約数・最小公倍数の間にはこのような関係式はありません）。
したがって、$\mathrm{lcm}(a,b)$ も $O(\log(a+b))$ で計算できることが分かります。
def lcm(a,b):
return a // gcd(a,b) * bgcdは当然 $a,b$ の約数ですから、積より前に除算を入れることが出来ます。特に64ビット整数などを扱っている場合には、先に積をとってしまうと（最小公倍数が64ビット整数の範囲内でとれる場合でも）桁あふれしてしまう可能性があります。
関連：ax+by=dの整数解
$d = \mathrm{gcd}(a,b)$ とするとき、$ax + by = d$ の解の1つをEuclidの互除法により計算できます。このことは高校数学でもしばしば取り上げられます。「拡張互除法」ということもあるようです。
def extgcd(a, b):
"""
ax + by = gcd(a,b) = d となる (x,y,d) を返す
"""
if b == 0:
return (1, 0, a)
q, r = a // b, a % b
x, y, d = extgcd(b, r)
s, t = y, x - q * y
return s, t, d 実装の際には、$s, t$ の式を正しく書けるかがポイントとなるでしょう。ここは
$bx + ry = d$ となる $x,y$ が与えられたとき、$(qb+r)s + bt = d$ となる $s,t$ を求めよ
という問題を考えればよいです。係数比較により $x = qs + t$ かつ $y = s$ が成り立つようにすればよいので、上で実装した通りとなります。なお、$x,y$ としては負の数が返される場合がある（一般的に、片方が正で片方が負）ので、用途によっては返り値の扱いに注意する必要があります。
関連：a mod Nの逆元
$a,N$ が整数で、$\mathrm{gcd}(a,N) = 1$ のとき、$ax \equiv 1\pmod{N}$ となる $x$ を $a\bmod N$ の逆元といいます。上で述べたように、Euclidの互除法から $ax + Ny = 1$ を満たす整数 $(x,y)$ が計算できます。この $x$ は当然 $ax\equiv 1\pmod{N}$ を満たすので、これが逆元です。$x < 0$ になって欲しくない場合は、負数ならば $N$ を加える処理を追加すればよいです。
def inv_mod(a, N):
"""
ax = 1 mod N となる x を返す
"""
return extgcd(a, N)[0]なお $N$ が素数の場合には、Fermatの小定理より、$a$の逆元を$a^{N-2} \bmod N$として計算できます。$a < N$ とすると、互除法による計算では計算量が$O(\log a)$、Fermatの小定理の方法では $O(\log N)$ となります。使い分ける必要はほとんどないと思いますが、小さい $a$ しか相手にする必要がない場合には互除法を使う方が計算量が少ないです（応用すると、単一の$a$ ではなく$1$ 以上 $n$ 以下のすべての $a$ の逆元を計算する必要がある場合に、$O(n)$ で計算することもできます）。
関連：中国剰余定理
$a_1$, $a_2$を整数として、$m_1$, $m_2$ を互いに素な整数とします。このとき
\[ x \equiv a_1\pmod{m_1},\qquad x\equiv a_2\pmod{m_2} \]を満たす整数 $x$ が、$\bmod m_1m_2$ で一意に定まります。このことを中国剰余定理といいます。
$x$ を計算してみましょう。第1の条件が成り立つような $x$ は、$x = a_1 + m_1y$ と書けます。上手く $y$ を調整して、$x\equiv a_2\pmod{m_2}$ が成り立つようにしましょう。
\[ x \equiv a_2\pmod{m_2} \iff a_1 + m_1y\equiv a_2\pmod{m_2} \iff m_1y\equiv a_2-a_1\pmod{m_2} \] となります。したがって、$m_1 \bmod m_2$ の逆元を $z$ とすれば、$y = (a_2 – a_1)z$ が条件を満たします。
def CRT(a_1, a_2, m_1, m_2):
"""
x = a_1 mod m_1, x = a_2 mod m_2 となる x を返す
"""
y = (a_2 - a_1) * inv_mod(m_1, m_2) % m_2
return a_1 + m_1 * y