中國剩余定理(解一元線性同余方程組)


問題:

  在《孫子算經》中有這樣一個問題:“今有物不知其數,三三數之剩二(除以3余2),五五數之剩三(除以5余3),七七數之剩二(除以7余2),問物幾何?”這個問題稱為“孫子問題”,該問題的一般解法國際上稱為“中國剩余定理”。

解析:

  題目意思即為有這樣一組方程:

(m1---mn兩兩互質)

 

 

 設:

 

 

 設:

 

 

 設:  所以:

 

 

 所以方程有解為:

 

 

 查找modM意義下,或最小正整數解:方程有唯一解:

注釋:假設整數m1,m2,m3...mn兩兩互質,則對於任意整數a1,a2,a3...an,x有解

例:

對於

x % 3 = 2

x % 5 = 3

x % 7 = 2

我們需要構造一個答案         //inv為求逆元

5*7*inv(5*7,  3) % 3  =  1

3*7*inv(3*7,  5) % 5  =  1

3*5*inv(3*5,  7) % 7  =  1 

然后兩邊同乘你需要的數,得

2 * 5*7*inv(5*7,  3) % 3  =  2

3 * 3*7*inv(3*7,  5) % 5  =  3

2 * 3*5*inv(3*5,  7) % 7  =  2

令 

a = 2 * 5*7*inv(5*7,  3) 

b = 3 * 3*7*inv(3*7,  5) 

c = 2 * 3*5*inv(3*5,  7) 

那么

a % 3 = 2

b % 5 = 3

c % 7 = 2

其實答案就是a+b+c

因為

a%5 = a%7 = 0 //a是5*7*inv(5*7,3)的倍數,所以模5和7都是

b%3 = b%7 = 0 //同理,b,c模3/7,3/5都是0

c%3 = c%5 = 0 

所以

(a + b + c) % 3 = (a % 3) + (b % 3) + (c % 3) = 2 + 0 + 0 = 2

(a + b + c) % 5 = (a % 5) + (b % 5) + (c % 5) = 0 + 3 + 0 = 3

(a + b + c) % 7 = (a % 7) + (b % 7) + (c % 7) = 0 + 0 + 2 = 2

即a+b+c是一個滿足條件得答案      //最小得正整數解為(a+b+c)%105(105=3*5*7,兩兩互質,105為最小公倍數)

代碼:

#include <algorithm>
#include <iostream>
#include <cstring>
#include <cstdlib>
#include <vector>
#include <cstdio>
#include <queue>
#include <cmath>
#include <map>
#include <set>

using namespace std;

typedef long long ll;

const int maxn=1e5+10;

ll a[maxn],m[maxn],n;

ll ex_gcd(ll a,ll b,ll &x,ll &y){
    if(b==0){x=1;y=0;return a;}
    ll ans=ex_gcd(b,a%b,x,y);
    ll temp=x;
    x=y;
    y=temp-a/b*y;
    return ans;
}
ll inv(ll a,ll b){   //求逆
   ll x,y;
   ll ans=ex_gcd(a,b,x,y);
   if(ans!=1)return -1;
   if(x<0)x=(x%b+b)%b;
   return x;
}
ll China(){//中國剩余定理
   ll M=1;
   for(int i = 0;i<n;i++){
      M*=m[i];
   }
   ll sum=0;
   for(int i=0;i<n;i++){
      ll res=M/m[i];
      sum=(sum+a[i]*res*inv(res,m[i]))%M;
   }
   return sum;
}

int main(){
    while(scanf("%lld",&n)!=EOF){
        for(int i=0;i<n;i++){
            scanf("%lld%lld",&m[i],&a[i]);
        }
        ll ans=China();
        printf("%lld\n",ans);
    }
    return 0;
}

上述方法只能解決m兩兩互質的同余方程組,如果不互質呢?

法一:

假設只有兩個方程,如果能用一個方程代替兩個方程,並能得到這個方程的解,那么這個問題就解決了:

x=a1+ k1*b1 (x%b1= a1)

x=a2+ k2*b2 (x%b2 = a2)

則: a1l+ k1*b1= a2 + k2*b2

then: b1* k1 + b2*(-k2)= a2-a1

假設gcd(b1,b2)= g若(a2 - a1)%g!=0 :無解

否則:k1*b1= (a2-a1) + k2*b2

then: k1*b1/g = (a2-a1)/g + k2*b2/g

then: k1*b1/g = (a2-a1)/g (mod b2/g)

then : k1 = inv(b1/g,b2/g)*(a2-a1)/g (mod b2/g)

then:x = a1 + inv(b1/g,b2/g)*(a2-a1)/g*b1 + b2*b1/g*y

then :x = c(mod m)

中:c= a1 + inv(b1/g,b2/g)*(a2-a1)/g*b1

m= b1*b2/g = lcm(b1,b2)

當作一個新方程繼續與下一個聯立,直到剩下最后一個方程xn = cn(mod mn)則,

任意一個xn都是滿足條件的解,取最小正整數即可ans = (cn%mn + mn)%mn

代碼:

#include <algorithm>
#include <iostream>
#include <cstring>
#include <cstdlib>
#include <vector>
#include <cstdio>
#include <queue>
#include <cmath>
#include <map>
#include <set>

using namespace std;

typedef long long ll;
const int maxn=1e5+10;
ll C[maxn],M[maxn];
ll n,k;
ll gcd(ll a,ll b){
    return b==0?a:gcd(b,a%b);
}
ll ex_gcd(ll a,ll b,ll &x,ll &y){
    if(b==0){x=1;y=0;return a;}
    ll g=ex_gcd(b,a%b,x,y);
    ll temp=x;
    x=y;
    y=temp-a/b*y;
    return g;
}
ll inv(ll a,ll mod){
    ll X,Y;
    ll g=ex_gcd(a,mod,X,Y);
    if(g!=1)return -1;
    return (X%mod+mod)%mod;
}
/*
ll mul(ll a,ll b,ll mod){   //快速乘法
    ll ans=0;
    while(b){
        if(b&1)ans=(ans%mod+a%mod)%mod;
        b>>=1;
        a=(a%mod+a%mod)%mod;
    }
    return ans;
}
*/
int main(){
    while(scanf("%lld",&n)!=EOF){
        for(int i = 0;i<n;i++){
            scanf("%lld%lld",&M[i],&C[i]);
        }
        bool flag=true;
        for(int i=1;i<n;i++){
            ll M1=M[i-1],M2=M[i],C1=C[i-1],C2=C[i];
            ll g=gcd(M1,M2);
            if((C2-C1)%g){flag = false;break;}
            M[i]=M1/g*M2;   //可能會爆
            ll INV=inv(M1/g,M2/g);
            if(INV==-1){flag=false;break;}
            C[i]= C1+(INV*((C2-C1)/g))%(M2/g)*M1;
            C[i]=(C[i]%M[i]+M[i])%M[i];
        }
        if(!flag)printf("-1\n");
        else printf("%lld\n",C[n-1]);
    }
    return 0;
}

法二:(僅限m比較小的時候/牛客2019夏多校1-d)

注意構造前i個方程的最小自然數解,假設前i-1個方程的最小自然數解為f(i-1)。lcm(p1,p2,……,pi-1)=y,令f(i)=f(i-1)然后讓f(i)每次加y,知道滿足第i個方程即可。可以證明在保證有解的情況下,加法運算不會超過pi次。(所以要先判斷方程組有無解)

#include <algorithm>#include <iostream#include <cstring#include <cstdlib>#include <vector#include <cstdio>

#include <queue> #include <cmath> #include <map> #include <set> using namespace std; typedef long long ll; const int maxn=100+10; ll ans,n,a[maxn],p[maxn],base,M; const char *definitely_lie="he was definitely lying"; const char *probably_lie="he was probably lying"; int main(){ scanf("%lld%lld",&n,&M); for(int i=0;i<n;i++)scanf("%lld%lld",&p[i],&a[i]); for(int i=0;i<n;i++){ for(int j=i+1;j<n;j++){ ll d=__gcd(p[i],p[j]); if(d!=1){ if(a[i]%d!=a[j]%d){ puts(definitely_lie); return 0; } } } } ans=0; base=1; for(int i=0;i<n;i++){ while(ans%p[i]!=a[i]){ ans+=base; if(ans>M){ puts(probably_lie); return 0; } } ll gg=p[i]/__gcd(base,p[i]); if(M/base>=gg)base*=gg; else{ for(int j=i+1;j<n;j++)if(ans%p[j]!=a[j]){ puts(probably_lie); return 0; } printf("%lld\n",ans); return 0; } }
   for(int i=0;i<n;i++)assert(ans%p[i]==a[i]);
   printf("%lld\n",ans);
return 0; }

 


免責聲明!

本站轉載的文章為個人學習借鑒使用,本站對版權不負任何法律責任。如果侵犯了您的隱私權益,請聯系本站郵箱yoyou2525@163.com刪除。



 
粵ICP備18138465號   © 2018-2025 CODEPRJ.COM