中国剩余定理(CRT)：求解模线性方程组

定义

 设有同余式组 (S)：

(S)：⎧⎩⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪x≡a1(modm1)x≡a2(modm2)x≡a3(modm3)⋮x≡ak(modmk)

即

x≡ai(modmi),i=1,2,..,k

 其中 ai 和 mi 已知，且当 i≠j 时， 有 (mi,mj)=1 , 求x对模 m1m2...mk 的解。

解的情况

 有唯一解， 证明略。

 下面会构造出此解，然后证明其正确性和唯一性。

 设

M=∏ki=1mi//即所有m的积

Mi=Mmi//即除mi以外所有m的积

//M−1i为Mi的逆元(对模mi)，即MiM−1i≡1(modmi)

 则(S)的解为:

x=tM+∑i=1kaiMiM−1i

即

x=(∑i=1kaiMiM−1i)(modM)

 下证其正确性:

我们考查任一在[1, k]内的第i项: aiMiM−1i , 有

aiMiM−1i≡{ai⋅1≡(modmj),j=i0≡(modmj),j≠i

这里不难理解， 第二个为0是因为 mj|Mi 。

所以x满足如下条件:

∀i∈{1,2,..,k} , 有

x=aiMiM−1i+∑j≠i(ajMjM−1J)≡ai+∑j≠i0≡ai(modmi)

综上， 此解正确。

 下面证明其唯一性

假设 x1 和 x2 都是方程组(S)的解，那么：

∀i∈{1,2,..,k} , 有 x1−x2≡0(modmi)

即 mi|(x1−x2)

那么显然 M|(x1−x2)

即 x1与x2 之间相差M的整数倍，那么在模M的范围内，只有一解。

更通用的解法

 前面说的模线性方程组是k个m互素的情况， 那么更一般的呢？

 更一般的模线性方程也能求解，只不过需要从前往后遍历每一个式子， 每一个都与下一个式子合并，作为一个新的式子继续与下一个合并， 最后得出结果。

 重点！！！

 我们先假设k=2时的情况，然后再慢慢推广。

 此时(S)为

(S)：{x≡a1(modm1)(1)x≡a2(modm2)(2)

即

{x=k1m1+a1(3)x=k2m2+a2(4)

 我们将(3)(4)联立，得

k1m1=k2m2+a2−a1

, 即

k1m1≡a2−a1(modm2)(5)

//这个式子特别重要

 因为(5)式未知数只有 k1 ，所以我们用扩展gcd可以求出 k1 的值， 不妨设求出的 k1=t , 则对某t’, 有 k1=t+t′m2 。我们将它代入(3)式，得到

x≡(a1+tm1)(modlcm(m1m2))(6)//lcm为最小公倍数

 那么(6)式中的x即为同时满足(1)和(2)的x值，至此k=2的情况讨论完毕。

 至于k>2的情况， 我们可以继续将(6)式与(3)式像(1)式和(2)式那样执行上面的操作，这样到最后得到的x就满足所有的条件了。

 具体的操作可参考代码理解， 这里简单说一下中间量的维护：

 每次对2个式子进行以上操作时， 第二个式子和以前一样，不作讨论，而第一个式了如果要写成 x≡a1(modm1) 的形式的话， a1 和 m1 都是要变的，那么变成什么样子决定着上述算法的最终实现。

 首先来说 m1 的变化 ， 从(1)到(6)， m1 变成了lcm( m1,m2 ), 如果接着往下走的话， m1 又会变成lcm( m1,m2,m3 ),…, 直到最后， m1变成所有m的最小公倍数。

 而对于 a1 的变化 ，我们可以这样想，如果k=1的话， 结果ans是等于 a1 的，而k = 2的时候， 观察(6)我们会发现， a1 变成了 a1+tm1 , 那结果呢， 结果ans(不看模)也变成了 a1+tm1 , 如果再往下写你会发现， ans 不是变成了 a1+tm1 ， 而是变成了 ans+t⋅lcm , 其中t 是由我们解(5)得到的。

 这里说的有点乱，具体的细节可以结合下面的实现代码去理解。

实现

 m[i] 与m[j] 互质:

int CRT(int a[], int m[], int k) {
    int i, d, x, y, Mi, ans = , M = ;
    for (i = ; i < k; i++) M *= m[i]; 
    for (i = ; i < k; i++) {
        Mi = M / m[i];
        d = extgcd(m[i], Mi, x, y);     // y 为逆元 -- Mi*y === 1 (% m[i])
        ans = (ans + a[i]*y*Mi) % M;
    }
    if (ans > ) return ans;
    else return (ans + M);
}
           

 m[i] 与m[j] 不一定互质:

LL mod(LL a, LL m) {    //处理取模
    LL res = a%m;
    if(res <= ) res += m;
    return res;
}
bool CRT(LL a[], LL m[], LL k, LL& ans) {
    ans = a[];
    LL lcm = m[], x, y, d;
    if(ans == ) ans = m[];
    for(LL i = ; i < k; ++i) {
        d = extgcd(lcm, m[i], x, y);    //求t: t = (a[i]-ans)/d*x;
        if((a[i]-ans)%d) return ;
        ans = mod(ans + lcm*mod((a[i]-ans)/d*x, m[i]/d),(lcm/d*m[i]));
        lcm = lcm/d*m[i];
    }
    return ;
}
           

例题

POJ - 1006 Biorhythms

HDU - 1573X问题

POJ - 2891 Strange Way to Express Integers

平生写过的最浮夸的博客。。