初学--求解模线性方程组(中国余数定理)。
2013-08-15 19:17
387 查看
中国剩余定理,孙子高富帅。
中国剩余定理到求解运用到了扩展欧几里得算法。
求解模线性方程组(中国余数定理)
a=B[1](mod W[1])
a=B[2](mod W[2])
........
a=B
(mod W
)
其中W,B已知,W[i]>0且W[i]与W[j]互质, 求a
设m1,m2,…mk是两两互素的正数,则对任意的整数b1,b2,…bk,同余方程组
x1 = b1 mod m1
x2 = b2 mod m2
…
xk = bk mod mk
其解为:
X = (M1’*M1*b1)+(M2’*M2*b2)+…+(Mk’*Mk*bk) mod m;
其中
m = m1*m2*…*mk;
Mi = m / mi;
Mi’是Mi的逆元
注意:这个求解的方法到前提是,w[i],w[j]是互质的。
具体到代码贴一下:(还是比较好理解到)
现在的问题是,当w[i] ,w[j] 不互质的时候,怎么去解决。这个也是我们要解决到问题了。
做到题目里,很少有全部互质到,那就很简单了,直接套模板了,这个你懂了,也是模板。
如何解决??????
合并操作.
先贴一段别人到思路:
假设C ≡ A1 (mod B1),C ≡ A2 (mod B2)。令C = A1 + X1B,那么X1B1 ≡ A2 − A1 (mod B2)。
用扩展欧几里德算法求出X1,也就求出C。令B = lcm(B1, B2),
那么上面两条方程就可以被C’ ≡ C (mod B)代替。迭代直到只剩下一条方程。
有点跳跃性:
C≡A1 (mod B1) ==>C=B1*k1 + A1; //这个容易理解的
C≡A2 (mod B2) ==>C=B2*k2 + A2;
合并一下 B1*k1 + A1 = B2*k2 + A2 ==> B1*k1-B2*k2 = A2-A1;
对比一下扩展欧几里得到式子 ax+by=c ;
则a=B1 ; b= -B2; c=A2-A1;
根据扩展欧几里得可以求的x0,y0;
gcd(B1,B2) = d;
c=A2-A1;
t=B2/d; //扩展欧几里得里也有这个 b=b/d;
x到最小非负解为 x=(x*c/d %t +t)%t; //需要理解一下。可以分两步 x=x*c/d; x=(x%t +t )%t;
那此时就能求出C 了 C=B1*K1 + A1 ==> C=B1*x + A1; //就是带入方程,求出C { C’ ≡ C (mod B) }
B=(B1*B2)/d; //最小公倍数到求法 b=(b1*b2)/gcd(b1,b2);
具体到看一看代码。
中国剩余定理到求解运用到了扩展欧几里得算法。
求解模线性方程组(中国余数定理)
a=B[1](mod W[1])
a=B[2](mod W[2])
........
a=B
(mod W
)
其中W,B已知,W[i]>0且W[i]与W[j]互质, 求a
设m1,m2,…mk是两两互素的正数,则对任意的整数b1,b2,…bk,同余方程组
x1 = b1 mod m1
x2 = b2 mod m2
…
xk = bk mod mk
其解为:
X = (M1’*M1*b1)+(M2’*M2*b2)+…+(Mk’*Mk*bk) mod m;
其中
m = m1*m2*…*mk;
Mi = m / mi;
Mi’是Mi的逆元
注意:这个求解的方法到前提是,w[i],w[j]是互质的。
具体到代码贴一下:(还是比较好理解到)
int china(int b[],int m[],int len) { int i,d,x,y,mi,sum=1,hxl=0;//好多变量 for(i=1;i<=len;i++) sum=sum*m[i];//sum=m1*m2*m3*.....mn for(i=1;i<=len;i++) { mi=sum/m[i];//Mi=m/mi d=Ex_gcd(m[i],mi,x,y);//y是mi的逆元,为什么是y呢,我感到很困惑。 hxl=(hxl+mi*y*b[i])%sum;//.. } if(hxl>0) return hxl; else return hxl+sum; }
现在的问题是,当w[i] ,w[j] 不互质的时候,怎么去解决。这个也是我们要解决到问题了。
做到题目里,很少有全部互质到,那就很简单了,直接套模板了,这个你懂了,也是模板。
如何解决??????
合并操作.
先贴一段别人到思路:
假设C ≡ A1 (mod B1),C ≡ A2 (mod B2)。令C = A1 + X1B,那么X1B1 ≡ A2 − A1 (mod B2)。
用扩展欧几里德算法求出X1,也就求出C。令B = lcm(B1, B2),
那么上面两条方程就可以被C’ ≡ C (mod B)代替。迭代直到只剩下一条方程。
有点跳跃性:
C≡A1 (mod B1) ==>C=B1*k1 + A1; //这个容易理解的
C≡A2 (mod B2) ==>C=B2*k2 + A2;
合并一下 B1*k1 + A1 = B2*k2 + A2 ==> B1*k1-B2*k2 = A2-A1;
对比一下扩展欧几里得到式子 ax+by=c ;
则a=B1 ; b= -B2; c=A2-A1;
根据扩展欧几里得可以求的x0,y0;
gcd(B1,B2) = d;
c=A2-A1;
t=B2/d; //扩展欧几里得里也有这个 b=b/d;
x到最小非负解为 x=(x*c/d %t +t)%t; //需要理解一下。可以分两步 x=x*c/d; x=(x%t +t )%t;
那此时就能求出C 了 C=B1*K1 + A1 ==> C=B1*x + A1; //就是带入方程,求出C { C’ ≡ C (mod B) }
B=(B1*B2)/d; //最小公倍数到求法 b=(b1*b2)/gcd(b1,b2);
具体到看一看代码。
#include <iostream> #include <cstdio> #include <cstdlib> #include <cstring> using namespace std; __int64 x, y, t; __int64 egcd(__int64 a, __int64 b) { if (b == 0) { x = 1; y = 0; return a; } else { __int64 e = egcd(b, a % b); t = x; x = y; y = t - a / b * y; return e; } } int main() { //freopen("in.txt", "r", stdin); __int64 m1, m2, r1, r2, d, c, t; bool sb; int n; while (cin >> n) { sb = 0; cin >> m1 >> r1; for (int i = 0; i < n - 1; i++) { cin >> m2 >> r2; if (sb) continue;//一旦有无解情况出现,自然就不用在计算了 d = egcd(m1, m2);//Ex_gcd(a,b); c = r2 - r1;//ax+by=(c2-c1) c=c2-c1;//c 是余数 if (c % d)//是否有解 { sb = 1; continue; } t = m2 / d;//b= b/d; x = (c / d * x % t + t) % t;//x=(x0*c/d % b+b)%b; //求最小到x r1 = m1 * x + r1;//C = A1 + X1B1 ==> A1=C; m1 = m1 * m2 / d;//最小公倍数赋值给m1 C’ ≡ C (mod B) } if (sb) cout << -1 << endl; else cout << r1 << endl; } }
相关文章推荐
- 中国剩余定理求解同余线性方程组(模数互素和非互素的情况)
- 中国剩余定理求解同余线性方程组—(互素和非互素的情况)
- 中国剩余定理(CRT):求解模线性方程组
- 中国剩余定理求解同余线性方程组—(互素和非互素的情况)
- Strange Way to Express Integers(扩展欧几里得+乘法逆元+中国剩余定理求解非互质的模线性方程组)
- HDU 1573 X问题 (中国剩余定理 模线性方程组)
- 数值分析 Gauss-Seidel迭代法求解线性方程组 MATLAB程序实现
- 算法导论-----数论-----中国余数定理
- 分别用雅可比(Jacobi)迭代法和高斯—塞德尔(Gauss—Seidel)迭代法求解线性方程组
- 初学A*算法求解静态地图的最短路径
- 采用GAUSS列主消元法求解线性方程组(MATLAB)
- 模线性方程组(中国剩余定理)
- 利用R和Octave求解线性方程组
- 利用中国剩余定理(求解一元线性同余方程组)
- lapack 线性方程组求解 C接口
- 高斯消元法求解线性方程组
- 数值计算-线性方程组求解(0)
- Problem 1402 猪的安家(中国余数定理)
- 中国剩余定理(余数定理)Chinese remainder
- java初学 中国象棋 总结