逆元の計算
合同式\(ax \equiv b \pmod{m}\)を解くことは、\(ax + my = b\)と等価です。これを拡張ユークリッドの互除法で解けます。
特に、\(a\)と\(p\)が互いに素の場合、\(a\)の\(p\)を法とする逆元を求めることができます。方法は\(ax \equiv 1 \pmod{p}\)を解くことです。
中国剰余定理（CRT）
問題設定
\(m_1, m_2, \ldots, m_n\)が\(n\)個の互いに素な整数であり、\(a_1, a_2, \ldots, a_n\)が与えられたとき、以下の連立合同式を解きます：
\( \begin{cases} x \equiv a_1 \pmod {m_1} \\ x \equiv a_2 \pmod {m_2} \\ \vdots \\ x \equiv a_n \pmod {m_n} \end{cases} \)
定理
\(M = \prod_{i=1}^n m_i\)（すべての\(m_i\)の積）、\(M_i = M/m_i\)（\(m_i\)以外のすべての\(m_j\)の積）、そして\(t_i\)を方程式\(M_i t_i \equiv 1 \pmod{m_i}\)の解（\(M_i\)の\(m_i\)を法とする逆元）とします。
このとき、解は\(x \equiv \sum_{i=1}^n a_i M_i t_i \pmod{M}\)で与えられます。
証明
\(i\)番目の合同式について、\(k \in [1,n]\)かつ\(k \neq i\)のとき、\(M_i\)の定義から\(a_i M_i t_i \equiv 0 \pmod{m_k}\)となります。また、\(t_i\)の定義から\(a_i M_i t_i \equiv a_i \pmod{m_i}\)となります。したがって、\(x = \sum_{i=1}^n a_i M_i t_i\)はすべての合同式を満たします。
一般解
整数\(k\)に対して、\(x + kM\)もまた解となります。
実装例：中国剰余定理の実装
#include <bits/stdc++.h>
using namespace std;
typedef long long i64;
i64 extended_gcd(i64 a, i64 b, i64 &x, i64 &y) {
if (b == 0) {
x = 1, y = 0;
return a;
}
i64 d = extended_gcd(b, a % b, x, y);
i64 temp = x - (a / b) * y;
x = y, y = temp;
return d;
}
int main() {
ios::sync_with_stdio(false);
cin.tie(nullptr);
int equations;
cin >> equations;
vector<i64> modulus(equations), remainder(equations);
i64 product = 1;
for (int i = 0; i < equations; i++) {
cin >> modulus[i] >> remainder[i];
product *= modulus[i];
}
i64 solution = 0;
for (int i = 0; i < equations; i++) {
i64 partial = product / modulus[i];
i64 inv, temp;
extended_gcd(partial, modulus[i], inv, temp);
inv = (inv % modulus[i] + modulus[i]) % modulus[i];
solution = (solution + remainder[i] * partial % product * inv % product + product) % product;
}
cout << solution << "\n";
return 0;
}
拡張中国剰余定理（EXCRT）
中国剰余定理は便利ですが、すべての法が互いに素である必要があります。より一般的なケースである法が互いに素でない場合を扱うために、拡張中国剰余定理が登場します。
問題
法が互いに素でない場合の連立合同式を解きます： \(\begin{cases} x \equiv r_1 \pmod {m_1} \\ x \equiv r_2 \pmod {m_2} \\ \vdots \\ x \equiv r_n \pmod {m_n} \end{cases}\)
導出
二つの合同式\(\begin{cases} x \equiv r_1 \pmod {m_1} \\ x \equiv r_2 \pmod {m_2} \end{cases}\)をマージすることを考えます。
これを等式の形で書くと： \(\begin{cases} x = r_1 + k_1 m_1 \\ x = r_2 + k_2 m_2 \end{cases}\)
連立して：\(r_1 + k_1 m_1 = r_2 + k_2 m_2\)
整理して：\(k_1 m_1 = (r_2 - r_1) + k_2 m_2\)
\(d = \gcd(m_1, m_2)\)とすると、両辺を\(d\)で割って： \(k_1 \frac{m_1}{d} = \frac{r_2 - r_1}{d} + k_2 \frac{m_2}{d}\)
\(\frac{m_2}{d}\)を法として考えると： \(k_1 \frac{m_1}{d} \equiv \frac{r_2 - r_1}{d} \pmod{\frac{m_2}{d}}\)
\(\frac{m_1}{d}\)と\(\frac{m_2}{d}\)は互いに素なので、逆元が存在します： \(k_1 \equiv \text{inv}\left(\frac{m_1}{d}, \frac{m_2}{d}\right) \cdot \frac{r_2 - r_1}{d} \pmod{\frac{m_2}{d}}\)
これを\(x = r_1 + k_1 m_1\)に代入して整理すると、新しい合同式が得られます。
注意すべき点：\(\frac{r_2 - r_1}{d}\)が整数である必要があります。つまり、\(d | (r_2 - r_1)\)でなければ解は存在しません。
実装
#include <bits/stdc++.h>
using namespace std;
typedef __int128 i128;
typedef long long i64;
i64 extended_gcd(i64 a, i64 b, i64 &x, i64 &y) {
if (b == 0) {
x = 1, y = 0;
return a;
}
i64 d = extended_gcd(b, a % b, x, y);
i64 temp = x - (a / b) * y;
x = y, y = temp;
return d;
}
i64 modular_inverse(i64 a, i64 mod) {
i64 x, y;
extended_gcd(a, mod, x, y);
return (x % mod + mod) % mod;
}
int main() {
ios::sync_with_stdio(false);
cin.tie(nullptr);
int num_equations;
cin >> num_equations;
vector<i64> mods(num_equations), rems(num_equations);
for (int i = 0; i < num_equations; i++) {
cin >> mods[i] >> rems[i];
}
i64 current_mod = mods[0];
i64 current_rem = rems[0];
for (int i = 1; i < num_equations; i++) {
i64 next_mod = mods[i];
i64 next_rem = rems[i];
i64 x, y;
i64 gcd_val = extended_gcd(current_mod, next_mod, x, y);
i64 diff = next_rem - current_rem;
if (diff % gcd_val != 0) {
cout << "-1\n";
return 0;
}
i64 lcm = current_mod / gcd_val * next_mod;
i64 temp = ((i128)diff / gcd_val * x % (next_mod / gcd_val) + (next_mod / gcd_val)) % (next_mod / gcd_val);
current_rem = current_mod * temp + current_rem;
current_rem = (current_rem % lcm + lcm) % lcm;
current_mod = lcm;
}
cout << current_rem << "\n";
return 0;
}