poj 2891 Strange Way to Express Integers (扩展欧几里得算法解一般线性同余方程组)

Posted by 111qqz on Friday, October 14, 2016

TOC

题目链接

题意:给出k个方程,形式为 x==r1,求最小的正数x,无解输出-1.

思路:首先很容易让人联想到crt.

然而crt的使用条件是,所有的m(也就是这道题中的a)两两互质,这道题并不满足,因此不能使用crt.

首先,我们看两个式子的情况
X mod m1=r1……………………………………………………………(1)
X mod m2=r2……………………………………………………………(2)
则有
X=m1k1+r1………………………………………………………………()
X=m2k2+r2
那么 m1
k1+r1=m2k2+r2
整理,得
m1
k1-m2*k2=r2-r1
令(a,b,x,y,m)=(m1,m2,k1,k2,r2-r1),原式变成
ax+by=m
熟悉吧?此时,因为GCD(a,b)=1不一定成立,GCD(a,b) | m 也就不一定成立。所以应该先判 若 GCD(a,b) | m 不成立,则方程无解。(理论依据:裴蜀定理
否则,继续往下。

解出(x,y),将k1=x反代回(),得到X。
于是X就是这两个方程的一个特解,通解就是 X'=X+k
LCM(m1,m2)
这个式子再一变形,得 X’ mod LCM(m1,m2)=X
这个方程一出来,说明我们实现了(1)(2)两个方程的合并。
令 M=LCM(m1,m2),R=r2-r1 注意这里原博客写错了,应该为R=x*m1+r1
就可将合并后的方程记为 X mod M = R。

然后,扩展到n个方程。
用合并后的方程再来和其他的方程按这样的方式进行合并,最后就能只剩下一个方程 X mod M=R,其中 M=LCM(m1,m2,…,mn)。
那么,X便是原模线性方程组的一个特解,通解为 X'=X+k*M。

如果,要得到X的最小正整数解,就还是原来那个方法:

X%=M;
if (X<0) X+=M;

参考博客

这篇题解是我看了那么多讲得最清楚的一篇了….虽然证明的那个地方笔误了。。。好在对照下代码就看出来了(个鬼,我问了好几个人,想到现在。。。最后发现是笔误。好气啊

这道题实际上给出了解线性同余方程组的一般做法(即,不要求所有的m两两互质)

这种做法的理论依据是扩展欧几里得算法,和中国剩余定理没什么关系…

不过既然都是解线性同余方程组。。。非要把这种做法叫什么ex_crt…我也不是很反对….

/* ***********************************************
Author :111qqz
Created Time :Sat 15 Oct 2016 02:42:54 AM CST
File Name :code/poj/2891.cpp
************************************************ */

#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
#include <vector>
#include <queue>
#include <set>
#include <map>
#include <string>
#include <cmath>
#include <cstdlib>
#include <ctime>
#define fst first
#define sec second
#define lson l,m,rt<<1
#define rson m+1,r,rt<<1|1
#define ms(a,x) memset(a,x,sizeof(a))
typedef long long LL;
#define pi pair < int ,int >
#define MP make_pair

using namespace std;
const double eps = 1E-8;
const int dx4[4]={1,0,0,-1};
const int dy4[4]={0,-1,1,0};
const int inf = 0x3f3f3f3f;
const int N=1E6+7;
int n;
LL a[N],r[N];

LL exgcd( LL a,LL b,LL &x,LL &y)
{
    if (b==0)
    {
    x = 1;
    y = 0;
    return a;
    }
    LL ret = exgcd(b,a%b,y,x);
    y-=x*(a/b); //简化版本的exgcd...
    return ret;
}
LL ex_crt( LL *m, LL *r,int n)
{
    LL M = m[1],R = r[1],x,y,gcd;
    for ( int i = 2 ; i <= n ; i++)
    {
    gcd = exgcd(M,m[i],x,y);
    if ((r[i]-R)%gcd) return -1;
    LL gx = m[i]/gcd;
    x = x*(r[i]-R)/gcd;
    x = x % gx;
    R = R + x*M;
    M = M / gcd * m[i];
    R%=M;
    }
    return R>0?R:R+M;
}

int main()
{
    #ifndef  ONLINE_JUDGE 
    freopen("code/in.txt","r",stdin);
  #endif

    while (~scanf("%d",&n))
    {
        for ( int i = 1 ; i <= n ; i++) scanf("%lld %lld",&a[i],&r[i]);
        printf("%lld\n",ex_crt(a,r,n));
    }

  #ifndef ONLINE_JUDGE  
  fclose(stdin);
  #endif
    return 0;
}

「真诚赞赏,手留余香」

111qqz的小窝

真诚赞赏,手留余香

使用微信扫描二维码完成支付