0. 介绍
主要通过混合积的表示来逐步求得同余方程的解。
对于同余方程
{
x
≡
v
0
(
m
o
d
m
0
)
x
≡
v
1
(
m
o
d
m
1
)
⋯
x
≡
v
k
−
1
(
m
o
d
m
k
−
1
)
\begin{equation*} \begin{cases} x \equiv v_0 \quad (\ \bmod \ m_0)\\ x \equiv v_1 \quad (\ \bmod \ m_1)\\ \cdots\\ x \equiv v_{k-1} \quad (\ \bmod \ m_{k-1})\\ \end{cases} \end{equation*}
⎩
⎨
⎧x≡v0( mod m0)x≡v1( mod m1)⋯x≡vk−1( mod mk−1)
满足
gcd
(
m
0
,
m
1
,
⋯
,
m
k
−
1
)
=
1
\gcd(m_0,m_1, \cdots, m_{k-1})=1
gcd(m0,m1,⋯,mk−1)=1
1. 传统构造方法
M = ∏ i = 0 k − 1 m i b i = M m i , b i × i n v i ≡ 1 ( m o d m i ) x = ∑ i = 0 k − 1 v i × b i × i n v i M=\prod_{i=0}^{k-1} m_i\\ b_i=\frac{M}{m_i}, b_i\times inv_i \equiv 1 \quad(\ \bmod\ m_i)\\ x=\sum_{i=0}^{k-1} v_i \times b_i \times inv_i M=i=0∏k−1mibi=miM,bi×invi≡1( mod mi)x=i=0∑k−1vi×bi×invi
这样构造保证了 x x x的存在性,但是在计算时 b i b_i bi比较大。
下面是示例代码
template<typename T>
T exgcd(T a, T b, T &x, T &y) {
if ( b == 0) {
x = 1;
y = 0;
return a;
}
T x_0;
T y_0;
T d = exgcd( b, a % b, x_0, y_0);
x = y_0;
y = x_0 - (a / b) * y_0;
return d;
}
ll CRT(int cnt, ll *a, ll *b) {
ll pia = 1;
ll res = 0;
for (int i = 0;i < cnt; i++) {
pia *= a[i];
}
for (int i = 0;i < cnt; i++) {
ll m = pia / a[i];
ll inv, tmp;
exgcd( m, a[i], inv, tmp);
inv = (inv % a[i] + a[i])% a[i];
__int128 sb = m * inv % pia * b[i] % pia;
sb %= pia;
res = (res + (ll) sb ) % pia;
}
return res;
}
2. Garner‘s算法
由中国剩余定理,我们知道必然存在 x < ∏ i = 0 k − 1 m i x <\prod_{i=0}^{k-1} m_i x<∏i=0k−1mi满足上面的同余方程组。
x
x
x的混合积表示为
x
=
t
0
+
t
1
m
0
+
⋯
+
t
k
−
1
∏
i
=
0
k
−
2
m
i
x=t_0 + t_1m_0+ \cdots+t_{k-1} \prod_{i=0}^{k-2} m_i
x=t0+t1m0+⋯+tk−1i=0∏k−2mi
不断对
x
x
x和
∏
i
=
0
s
m
i
\prod_{i=0}^{s} m_i
∏i=0smi应用带余除数法可以证明
t
i
t_i
ti的唯一性。
代入第一同余方程得到
t
0
≡
v
0
(
m
o
d
m
0
)
t_0\equiv v_0 \quad (\ \bmod \ m_0)
t0≡v0( mod m0)
取
t
0
=
v
0
t_0 =v_0
t0=v0, 将
x
x
x代入第二个方程
t
0
+
t
1
m
0
≡
v
1
(
m
o
d
m
1
)
t
1
m
0
≡
v
1
−
t
0
(
m
o
d
m
1
)
t
1
≡
(
v
1
−
t
0
)
m
0
−
1
(
m
o
d
m
1
)
t_0 + t_1m_0 \equiv v_1 \quad (\ \bmod \ m_1)\\ t_1m_0\equiv v_1-t_0 \quad (\ \bmod \ m_1) \\ t_1 \equiv(v_1-t_0) m_0^{-1} \quad (\ \bmod\ m_1)
t0+t1m0≡v1( mod m1)t1m0≡v1−t0( mod m1)t1≡(v1−t0)m0−1( mod m1)
我们令
X
j
=
∑
s
=
0
j
−
1
t
s
∏
a
=
0
s
−
1
m
a
X_j=\sum_{s=0}^{j-1}t_s \prod_{a=0}^{s-1} m_a
Xj=∑s=0j−1ts∏a=0s−1ma,也就是
x
x
x的混合积表示中的前
j
j
j
项;此时求 x x x相当于,求 X k X_k Xk。
在求 X j X_j Xj时,我们必然已经求得了 X j − 1 X_{j-1} Xj−1。
代入到第 j j j个同余方程,可以得到
X j − 1 + t j ( ∏ s = 0 j − 1 m s ) ≡ v j ( m o d m j ) t j ( ∏ s = 0 j − 1 m s ) ≡ ( v j − X j − 1 ) ( m o d m j ) t j ≡ ( v j − X j − 1 ) ( ∏ s = 0 j − 1 m s ) − 1 ( m o d m j ) X_{j-1} + t_j(\prod_{s=0}^{j-1}m_s) \equiv v_j \quad (\ \bmod \ m_j)\\ t_j(\prod_{s=0}^{j-1}m_s) \equiv (v_j - X_{j-1} ) \quad (\ \bmod \ m_j)\\ t_j \equiv(v_j-X_{j-1})(\prod_{s=0}^{j-1}m_s)^{-1} \quad (\ \bmod \ m_j) Xj−1+tj(s=0∏j−1ms)≡vj( mod mj)tj(s=0∏j−1ms)≡(vj−Xj−1)( mod mj)tj≡(vj−Xj−1)(s=0∏j−1ms)−1( mod mj)
因此我们可以递推的来求
x
x
x。
在应用密码学手册中给出了算法描述。
主要的解释上面已经有了,唯一需要注意的是求 ( ∏ i = 0 j − 1 m i ) − 1 ( m o d m j ) (\prod_{i=0}^{j-1} m_i) ^{-1} \quad(\ \bmod\ m_j) (∏i=0j−1mi)−1( mod mj)。
运用了乘积的逆等于逆的乘积这一性质, 证明如下
( ∏ i = 0 j − 1 m i ) ∏ i = 0 j − 1 ( m i ) − 1 ≡ ∏ i = 0 j − 1 m i m i − 1 ≡ 1 ( m o d m j ) (\prod_{i=0}^{j-1} m_i)\prod_{i=0}^{j-1}(m_i)^{-1} \equiv \prod_{i=0}^{j-1} m_i m_i^{-1} \equiv 1 \quad (\ \bmod\ m_j) (i=0∏j−1mi)i=0∏j−1(mi)−1≡i=0∏j−1mimi−1≡1( mod mj)
代码
#include <iostream>
#include <vector>
using ll = long long;
static constexpr int MAXN = 1e5+10;
ll a[MAXN];
ll b[MAXN];
template<typename T>
T exgcd(T a, T b, T &x, T &y) {
if ( b == 0) {
x = 1;
y = 0;
return a;
}
T x_0;
T y_0;
T d = exgcd( b, a % b, x_0, y_0);
x = y_0;
y = x_0 - (a / b) * y_0;
return d;
}
ll CRT_GARNER(int cnt, ll *m, ll *v)
{
std::vector<ll> C(cnt, 1);
for (int i = 1; i < cnt; i++) {
for (int j = 0; j < i; ++j) {
ll inv;
ll tmp;
exgcd(m[j], m[i], inv, tmp);
// printf("inv of m[%d] mod m[%d]: %d\n", j, i, inv);
while ( inv < 0)
inv += m[i];
C[ i ] = C[ i ] * inv % m[i];
}
}
ll u = v[ 0 ];
ll x = u;
ll M = m[0];
for (int i = 1;i < cnt; i++) {
u = (v[i] - x) * C[i] % m[i];
while (u < 0)
u += m[ i ];
x = x + u * M;
M *= m[i];
}
return x;
}
举个例子
{
x
≡
2
(
m
o
d
3
)
x
≡
3
(
m
o
d
5
)
x
≡
2
(
m
o
d
7
)
\begin{equation*} \begin{cases} x \equiv 2\quad(\ \bmod \ 3)\\ x \equiv 3 \quad (\ \bmod \ 5)\\ x \equiv 2 \quad (\ \bmod \ 7) \end{cases} \end{equation*}
⎩
⎨
⎧x≡2( mod 3)x≡3( mod 5)x≡2( mod 7)
将
x
x
x表示为混积的形式
x
=
t
0
+
2
t
1
+
15
t
2
x=t_0+2t_1+15t_2
x=t0+2t1+15t2
显然 t 0 = 2 t_0=2 t0=2;
代入到第二个同余方程
2
+
3
t
1
≡
3
(
m
o
d
5
)
3
t
1
≡
1
(
m
o
d
5
)
t
1
≡
3
−
1
(
m
o
d
5
)
2+3t_1 \equiv 3 \quad (\ \bmod \ 5)\\ 3t_1 \equiv1 \quad (\ \bmod \ 5)\\ t_1\equiv3^{-1} \quad (\ \bmod\ 5)\\
2+3t1≡3( mod 5)3t1≡1( mod 5)t1≡3−1( mod 5)
3
−
1
(
m
o
d
5
)
3^{-1} \quad (\ \bmod \ 5)
3−1( mod 5)可由扩展欧几里得算法求得
3
×
2
−
1
×
5
=
1
3
−
1
(
m
o
d
5
)
=
2
3 \times 2 -1 \times 5=1\\ 3^{-1} \quad( \ \mod \ 5) =2
3×2−1×5=13−1( mod 5)=2
因此
t
2
=
2
t_2=2
t2=2。
代入第三个同余方程
2
+
3
×
2
+
15
t
2
=
8
+
15
t
2
15
t
2
+
8
≡
2
(
m
o
d
7
)
15
t
2
≡
1
(
m
o
d
7
)
t
2
≡
15
−
1
(
m
o
d
7
)
2+3\times2+15t_2=8+15t_2\\ 15t_2+8 \equiv2 \quad (\ \bmod\ 7)\\ 15t_2 \equiv1 \quad ( \bmod\ 7)\\ t_2\equiv15^{-1}\quad(\ \bmod \ 7)
2+3×2+15t2=8+15t215t2+8≡2( mod 7)15t2≡1(mod 7)t2≡15−1( mod 7)
同样是利用扩展欧几里得求出逆元
t
2
≡
1
(
m
o
d
7
)
t_2 \equiv 1 \quad (\ \bmod \ 7)
t2≡1( mod 7)
代入到混合积形式
x
=
8
+
15
=
23
x=8+ 15=23
x=8+15=23
参考
oi-wiki-crt
handbook-of-applied-crypho