模线性方程组

　　模线性方程组：x=a1(Mod m1)

                               x=a2(Mod m2)

                               ......

                               x=an(Mod mn)

　　这个问题的源自《孙子算经》，其中有这样一道算术题：“今有物不知其数，三三数之剩二，五五数之剩三，七七数之剩二，问物几何？”按照今天的话来说：一个数除以3余2，除以5余3，除以7余2，求这个数.这样的问题，也有人称为“韩信点兵”.它形成了一类问题，也就是初等数论中解同余式.这类问题的有解条件和解的方法被称为“中国剩余定理”，这是由中国人首先提出的.首先来看一个简单的例子：

① 有一个数，除以3余2，除以4余1，问这个数除以12余几？
解：除以3余2的数有：
　　2， 5， 8， 11，14， 17， 20， 23….
它们除以12的余数是：
　　2，5，8，11，2，5，8，11，….
除以4余1的数有：
　　1， 5， 9， 13， 17， 21， 25， 29，….
它们除以12的余数是：
　　1， 5， 9， 1， 5， 9，….
　　一个数除以12的余数是唯一的.上面两行余数中，只有5是共同的，因此这个数除以12的余数是5.
　　如果我们把①的问题改变一下，不求被12除的余数，而是求这个数.很明显，满足条件的数是很多的，它是 5＋12×整数，
　　整数可以取0，1，2，…，无穷无尽.事实上，我们首先找出5后，注意到12是3与4的最小公倍数，再加上12的整数倍，就都是满足条件的数.这样就是把“除以3余2，除以4余1”两个条件合并成“除以12余5”一个条件.《孙子算经》提出的问题有三个条件，我们可以先把两个条件合并成一个.然后再与第三个条件合并，就可找到答案.
②一个数除以3余2，除以5余3，除以7余2，求符合条件的最小数.
解：　　当某数被3除余1对，即写上70（因为70是5和7的倍数，是3的倍数多1），余2时即写70×2＝140，这140仍是5和7的倍数，是3的倍数余2。某数被5除余1，即写上21（因为21是3和7的倍数、5的倍数余1），余2时，则写上21×2＝42，余3时，则写上21×3＝63。某数被7除余1时，即写上15（因为15是3和5的倍数，是7的倍数余1），余2时，则写上15×2＝30。根据题意，把70×2+21×2＋15×2计算出来结果。然后减去3、5、7的最小公倍数105，一直减到少于105为止，就得到了符合题目的数：
　　70×2＋21×3+15×2-105×2＝23
　　即此数是23。

　　通过这个例子，那么就容易退出互素的模方程的一般解法了：令M=m1*m2*...*mn，wi=M/mi，则gcd(wi,mi)=1，那么有wi*p=1(Mod mi) （其中p为未知数），则有 wi*p-mi*q=1，可以用扩展欧几里得算法求出p和q，设ei=wi*p，那么x=Σei*ai%M。

　　代码如下：

View Code

 1 LL a[N],m[N],phi[N];
 2 int n;
 3 
 4 LL china()
 5 {
 6     int i;
 7     LL M=1,w,d,y,x=0;
 8     for(i=0;i<n;i++)M*=m[i];
 9     for(i=0;i<n;i++){
10         w=M/m[i];
11         exgcd(m[i],w,d,d,y);
12         x=(x+y*a[i])%M;
13     }
14     return (x+M)%M;
15 }

　　

　　如果不互素怎么求呢？当然也可以模仿中国剩余定理的方法，两个方程依次求解为一个方程，也就是最开始的引例：除以3余2，除以4余1，转换成12的余数是5，具体的方法如下：

　　　x=a1(Mod m1)

           x=a2(Mod m2)

     则有 m1*x+a1=m2*y+a2  =>  m1*x-m2*y=a2-a1  从这里可以看出，如果方程要有解, d|(a2-a1)必须成立，这是模线性方程有解的充分必要条件。然后可以用扩展欧几里得求出x和y了，容易得出x=x*(a2-a1)/d，M=m1/gcd(m1,m2)*m2，这样一直合并递推下去就可了。

　　代码如下：

View Code

 1 LL a[N],m[N],phi[N];
 2 int n;
 3 
 4 LL Modline(int n)
 5 {
 6     LL d,x,y,A,M,Mod;
 7     A=a[n-1],M=m[n-1];
 8     n--;
 9     // m1*x-m2*y=a2-a1
10     while(n--){
11         exgcd(M,m[n],d,x,y);
12         if((A-a[n])%d!=0){
13             return -1;
14         }
15         Mod=m[n]/d;
16         x=(x*((a[n]-A)/d)%Mod+Mod)%Mod;
17         A+=M*x;
18         M=M/d*m[n];
19     }
20     return A;
21 }

 

转载于:https://www.cnblogs.com/zhsl/archive/2013/04/12/3017109.html