中国人剰余定理とは整数の剰余に関する定理の一つです。 ここでは、連立合同方程式を中国人剰余定理を使ったPythonプログラムで解く方法について説明します。

• 中国人剰余定理

  $k$ 個の整数 $n_1, n_2, …, n_k$ が互いに素ならば、任意の整数 $a_1, a_2, …, a_k$ に対して、次の連立合同方程式を満たす整数 $x$ が $n_1 n_2 \cdots{} n_k$ を法として一意的に存在する。

  \[\begin{aligned} x &\equiv a_1 \pmod{n_1} \\ x &\equiv a_2 \pmod{n_2} \\ &\;\;\vdots \\ x &\equiv a_k \pmod{n_k} \end{aligned}\]

  唯一の解 $x$ は、$1 \le i \le k$ について \(b_i = N / n_i, \; b_i' = b_i^{-1} \pmod{n_i}\) とおくと、次の式で求めることができる。

  \[x = \sum_{i=1}^k a_i b_i b'_i \pmod{N}\]

中国人剰余定理で連立合同方程式 $x \equiv a_i \pmod{n_i} \;\;\text{for}\; i = 1,…,k$ を解くためのアルゴリズムは以下の通りです。

• アルゴリズム

  1. $N = n_1 n_2 \cdots{} n_k$ を求める。
  2. \(i \in \{ 1,...,k \}\) に対して、整数 $n_i$ と $b_i = N / n_i$ を求める。この2つの整数は互いに素となる。
  3. $n_i$ を法とする $b_i$ の逆元 $b'_i = b_i^{-1} \pmod{n_i}$ を求める。

  4. 次の式を計算して解 $x$ を求める。

     \[x = \sum_{i=1}^k a_i b_i b'_i\]

  5. 最小の解 $x$ を求める。

     \[x \pmod{N}\]

このアルゴリズムをPythonで書くと次のようになります。

from functools import reduce

def xgcd(a, b):
    x0, y0, x1, y1 = 1, 0, 0, 1
    while b != 0:
        q, a, b = a // b, b, a % b
        x0, x1 = x1, x0 - q * x1
        y0, y1 = y1, y0 - q * y1
    return a, x0, y0

def modinv(a, m):
    g, x, y = xgcd(a, m)
    if g != 1:
        raise Exception('modular inverse does not exist')
    else:
        return x % m

def chinese_remainder(a, n):
    # a := [a1, a2, ..., ak]
    # n := [n1, n2, ..., nk]
    total = 0
    prod = reduce(lambda x, y: x*y, n)
    for n_i, a_i in zip(n, a):
        b_i = prod // n_i
        total += a_i * b_i * modinv(b_i, n_i)
    return total % prod

連立合同方程式を中国人剰余定理で解く問題として『孫子算経』の問題を使います。 つまり、次の連立合同方程式を解いてみます。

先ほどのプログラムを使って計算してみると、解は 23 となります。

# 連立合同方程式を中国人剰余定理で解く
a = [2, 3, 2]
n = [3, 5, 7]
print(chinese_remainder(a, n))
# => 23

『孫子算経』の問題の答えと一致しました。

以上です。

参考文献

• Number Theory - The Chinese Remainder Theorem
• Chinese remainder theorem - Rosetta Code
• 中国の剰余定理 - Wikipedia