模线性方程组

#1303 : 数论六·模线性方程组

时间限制: 10000ms

单点时限: 1000ms

内存限制: 256MB

描述

小Ho：今天我听到一个挺有意思的故事！

小Hi：什么故事啊？

小Ho：说秦末，刘邦的将军韩信带领1500名士兵经历了一场战斗，战死四百余人。韩信为了清点人数让士兵站成三人一排，多出来两人；站成五人一排，多出来四人；站成七人一排，多出来六人。韩信立刻就知道了剩余人数为1049人。

小Hi：韩信点兵嘛，这个故事很有名的。

小Ho：我觉得这里面一定有什么巧妙的计算方法！不然韩信不可能这么快计算出来。

小Hi：那我们不妨将这个故事的数学模型提取出来看看？

小Ho：好！

<小Ho稍微思考了一下>

小Ho：韩信是为了计算的是士兵的人数，那么我们设这个人数为x。三人成排，五人成排，七人成排，即x mod 3, x mod 5, x mod 7。也就是说我们可以列出一组方程：

x mod 3 = 2
x mod 5 = 4
x mod 7 = 6

韩信就是根据这个方程组，解出了x的值。

小Hi：嗯，就是这样！我们将这个方程组推广到一般形式：给定了n组除数m[i]和余数r[i]，通过这n组(m[i],r[i])求解一个x，使得x mod m[i] = r[i]。

小Ho：我怎么感觉这个方程组有固定的解法？

小Hi：这个方程组被称为模线性方程组。它确实有固定的解决方法。不过在我告诉你解法之前，你不如先自己想想怎么求解如何？

小Ho：好啊，让我先试试啊！

提示：模线性方程组

输入

第1行：1个正整数, N，2≤N≤1,000。

第2..N+1行：2个正整数, 第i+1行表示第i组m,r，2≤m≤20,000,000，0≤r<m。

计算过程中尽量使用64位整型。

输出

第1行：1个整数，表示满足要求的最小X，若无解输出-1。答案范围在64位整型内。

样例输入

3
3 2
5 3
7 2

样例输出

23

小Hi：一开始就直接求解多个方程不是太容易，我们从n=2开始递推：
已知：
x mod m[1] = r[1]
x mod m[2] = r[2]
根据这两个式子，我们存在两个整数k[1],k[2]：
x = m[1] * k[1] + r[1]
x = m[2] * k[2] + r[2]
由于两个值相等，因此我们有：
	m[1] * k[1] + r[1] = m[2] * k[2] + r[2]
=>	m[1] * k[1] - m[2] * k[2] = r[2] - r[1]
由于m[1],m[2],r[1],r[2]都是常数，若令A=m[1],B=m[2],C=r[2]-r[1],x=k[1],y=k[2]，则上式变为：Ax + By = C。
我们可以先通过gcd(m[1], m[2])能否整除r[2]-r[1]来判定是否存在解。
假设存在解，则我们通过扩展欧几里德求解出k[1],k[2]。
再把k[1]代入x = m[1] * k[1] + r[1]，就可以求解出x。
同时我们将这个x作为特解，可以扩展出一个解系：
X = x + k*lcm(m[1], m[2]) k为整数
lcm(a,b)表示a和b的最小公倍数。其求解公式为lcm(a,b)=a*b/gcd(a,b)。
将其改变形式为：
X mod lcm(m[1], m[2]) = x。
令M = lcm(m[1], m[2]), R = x，则有新的模方程X mod M = R。

#include <bits/stdc++.h>

using namespace std;
typedef long long LL;
struct node
{
    int m,r;
};
node num[1011];
LL m1,r1,m2,r2,x,y,d;
void exgcd(LL a,LL b)
{
    if(b==0)
    {
        x=1;
        y=0;
        d=a;
        return ;
    }
    exgcd(b,a%b);
    LL temp=x;
    x=y;
    y=temp-a/b*y;
    return ;
}
bool solve()
{
     exgcd(m1,m2);
     if((r2-r1)%d)return false;
     LL a=m1/d;
     LL b=m2/d;
     LL c=(r2-r1)/d;
     exgcd(a,b);
     x=((c*x)%b+b)%b;
     r1=x*m1+r1;
     m1=m1*b;
     r1=((r1%m1)+m1)%m1;
     return true;
}
int main()
{
    int n;
    scanf("%d",&n);
    for(int i=0;i<n;i++)
        scanf("%d%d",&num[i].m,&num[i].r);
    m1=num[0].m;
    r1=num[0].r;
    for(int i=1;i<n;i++)
    {
        m2=num[i].m;
        r2=num[i].r;
        if(!solve())
        {
            printf("-1\n");
            return 0;
        }
    }
    printf("%lld\n",r1);
    return 0;
}