数论 - 中国剩余定理（孙子定理）

一、起源

孙子定理是中国古代求解一次同余式组的方法。是数论中一个重要定理。又称中国余数定理或中国剩余定理。一元线性同余方程组问题最早可见于中国南北朝时期（公元5世纪）的数学著作《孙子算经》卷下第二十六题，叫做“物不知数”问题，原文如下：有物不知其数，三三数之剩二，五五数之剩三，七七数之剩二。问物几何？即，一个整数除以三余二，除以五余三，除以七余二，求这个整数。《孙子算经》中首次提到了同余方程组问题，以及以上具体问题的解法，因此在中文数学文献中也会将中国剩余定理称为孙子定理。

二、问题描述

给定 2n 个整数 a1,a2,…,an 和 m1,m2,…,m，求一个最小的非负整数 x，满足 ∀i∈[1,n],x≡mi(mod ai)。

输入格式

第 1 行包含整数 n。

第 2…n+1 行：每 i+1 行包含两个整数 ai 和 mi，数之间用空格隔开。

输出格式

输出最小非负整数 x，如果 x 不存在，则输出 −1。

数据范围

1≤ai≤2³¹−1,
0≤mi<ai
1≤n≤25
所有 mi 的最小公倍数在 64 位有符号整数范围内。

输入样例

2
8 7
11 9

输出样例

31

三、算法思路

• 我们需要求的是一个数x使得：x mod a₁ = m₁，x mod a₂ = m₂，x mod a₃ = m₃，……，x mod a_n = m_n
• 我们可以先取出两式子，看是否能转化为一个式子，取出x mod a₁ = m₁，x mod a₂ = m₂
• 这两个式子等价于：x = k₁a₁+m₁和x = k₂a₂+m₂
• 得到：k₁a₁+m₁ = k₂a₂+m₂
• 得到：k₁a₁ - k₂a₂ = m₂ - m₁ ，这就是上一篇博客《数论 - 欧拉函数、快速幂、扩展欧几里得算法》中所阐述的扩展欧几里得算法（如若不知，可先翻阅上一篇博客），此方程只有当(a₁,a₂)|(m₂-m₁)时有解，其中(a₁,a₂)表示两者的最大公约数
• 运用扩展欧几里得算法求出一组k₁和k₂后可得：所有的组解为k₁+ka₂/d和k₂+ka₁/d，k为任意整数，要注意还需要翻(m₂-m₁)/(a₁,a₂)倍
• 所以x = k₁a₁+m₁且k₁=k₁+ka₂/d，等价于x = a₁k₁ + m₁ + ka₁a₂/d
• 可以令k₁ = a₁k₁ + m₁ ，m₁ = a₁a₂/d，则两个式子x = k₁a₁+m₁和x = k₂a₂+m₂化为了一个式子x = k₁a₁+m₁，且结构完全一致，证明可以合并！！！
• 由此反复合并，最终所有式子都化为一个式子解出x的最小的非负整数解
• 注意：为了防止溢出，每次都要将k₁降到最小正整数解，同样最后的x也要进行此操作得最小的非负整数解。比如：要将a降为mod b的最小正整数，则 a = (a%b + a) %b。

四、代码

#include <iostream>

using namespace std;

typedef long long LL;

//扩展欧几里得算法求裴蜀定理系数x,y
LL exgcd(LL a, LL b, LL &x, LL &y)
{
    if (!b)
    {
        x = 1;
        y = 0;
        return a;
    }

    LL d = exgcd(b, a % b, y, x);
    y -= a / b * x;
    return d;
}

int main()
{
    int n;
    cin >> n;

    LL x = 0; //存储最终答案 如果没解则赋值为 -1
    LL a1, m1;
    cin >> a1 >> m1;

    for(int i = 0; i < n - 1; i ++)
    {
        LL a2, m2;
        cin >> a2 >> m2;

        LL k1, k2;
        LL d = exgcd(a1, a2, k1, k2);
        if((m2 - m1) % d) //判断是否无解
        {
            x = -1;
            break;
        }

        k1 *= (m2 - m1) / d; //翻倍
        k1 = (k1 % (a2 / d) + (a2 / d)) % (a2 / d); //保证k1为符合条件的最小值

        m1 = k1 * a1 + m1;
        a1 = abs(a1 / d * a2);
    }

    if(x != -1) x = ((m1 % a1) + a1) % a1; //保证x为符合条件的最小值

    cout << x << endl;

    return 0;
}