中国剩余定理(孙子定理)
“中国剩余定理”是公元5-6世纪、我国南北朝时期的一部著名算术著作《孙子算经》中的一个“物不知数”的解法问题:今有物不知其数,三三数之剩二,
五五数之剩三,七七数之剩二。问物几何?答曰:二十三。
1.除数两两互质的情况
根据上面我们可以得到一组式子:
X ≡ 2 (mod 3)
X ≡ 3 (mod 5)
X ≡ 2 (mod 7)
即:
X % 3 = 2
X % 5 = 3
X % 7 = 2
那么我们设 X = t1 * a + t2 * b + t3 * c
继续设
a ≡ 1 (mod 3) b ≡ 0 (mod 3) c ≡ 0 (mod 3)
a ≡ 0 (mod 5) b ≡ 1 (mod 5) c ≡ 0 (mod 5)
a ≡ 0 (mod 7) b ≡ 0 (mod 7) c ≡ 1 (mod 7)
a = LCM(5,7) * p(a%3=1) b = LCM(3,7)*p(b%5=1) c = LCM(3,5)*p(c%7=1)
= 35 * 2 = 70 = 21 * 1 = 21 = 15 * 1 = 15
那么我们就可以根据相对应的颜色得到t1 = 2, t2 = 3, t3 = 2
所以 X = 2a + 3b + 2c = 2 * 70 + 3 * 21 + 2 * 15 = 233
233是X的一个解,通解是 X = 233 + gcd(3,5,7)*k = 233 + 105*k
这些a、b、c(下面同称为Ni)该如何求呢(除数3、5、7统称为mi)
X ≡ y1 (mod m1)
X ≡ y2 (mod m2)
X ≡ y3 (mod m3)
...
X ≡ an (mod mn)
设X = N1 * y1 + N2 *y2 + N3 * y3 + ... + Nn * yn;
设 Ni能够被m1、m2、m3、.... 、mi-1、 mi + 1、... 、mn整除,但就是不能被mi整除且Ni%mi = 1
设 Ni = LCM(m1, m2, m3, ..., mn) / mi * x = mi * y + 1, 要求Ni就必须得求x
上面式子可以转化为 LCM(m1, m2, m3, ..., mn) / mi * x + (- mi) * y = 1
转化之后的式子就是扩展欧几里德了, 这样就可以求出x
模板:
int CRT(int a[], int b[], int n) {//a[i]表示除数,b[i]表示余数 int M = 1, N, ans = 0, x, y; for(int i = 0 ; i < n ; i++) M *= a[i];//因为除数两两互质,所以最小公倍数就是其乘积 for(int i = 0 ; i < n ; i++) { N = M / a[i]; gcd(N, a[i], x, y);//利用扩展欧几里德求x ans += N * x * b[i]; } return ans % M; }
2.除数不两两互质
当除数不两两互质时上面的做法就不正确了,那么应该怎们做呢
我们可以两个的合并方程,如下面两个方程:
X ≡ b1 (mod n1)
X ≡ b2 (mod n2)
X = b1 + n1 * k1; (1)
X = b2 + n2 * k2;
可将上面两个方程画等: b1 + n1 * k1 = b2 + n2 * k2 ==> n1 * k1 = b2 - b1 + n2 * k2 (2)
设r = gcd(n1, n2),左右两边同时除以r得:
b1/r + n1*k1/r = b2/r + n2*k2/r
转化为 n1*k1/r = b2/r - b1/r + n2*k2/r
所以 n1*k1/r = (b2 - b1)/r(mod n2/r)
==> n1*k1 = (b2 - b1)(mod n2/r)
k1 = (b2-b1)/n1(mod n2/r)
令k’ = (b2 - b1) / n1
所以 k1 = k’(mod n2/r)
k1 = k’ + n2/r * C
将k1代入(1)式中:
X = b1 + n1 * k1 = b1 + n1 * (k’ + n2/r * C ) = b1 + n1 * k’ + n1 * n2 / r * C
由(2)得:
n1 * k1 = (b2 - b1)(mod n2)
所以 n1 * k1 + n2*k2 = b2 - b1
可用扩展欧几里德定理来求k1
t = n2 / r
k’ = (k1 % t + t) % t
其实k与k'是等价的
详解看下面链接:
http://972169909-qq-com.iteye.com/blog/1266328
模板:
ll CRT2(ll n[], ll b[], ll m) { int f = 0; ll n1 = n[0], n2, b1 = b[0], b2, c, t, k, x, y; for(ll i = 1 ; i < m ; i++) { n2 = n[i]; b2 = b[i]; c = b2 - b1; gcd(n1, n2, x, y);//扩展欧几里德 if(c % r != 0)//无解 { f = 1; break; } k = c / r * x;//扩展欧几里德求得k t = n2 / r; k = (k % t + t) % t;//这个k就是上面的k’ b1 = b1 + n1 * k; n1 = n1 * t; } if(f == 1) return -1; return b1; }