中国剩余定理

$$\{\begin{array}{l}x\equiv a_1(mod{m}_1)\\ x\equiv a_2(mod{m}_2)\\ ⋮\\ x\equiv a_r(mod{m}_r)\end{array}$$

问题

对于以上同余方程组，求解改同余方程组的未知数 $x$ 。

中国剩余定理(CRT)

设 $m_1,m_2,...,m_r$ 为两两互质的正整数，则同余方程组有整数解。

1. 所有模数的积 $M=m_1{m}_2\cdots m_r$
2. 计算第 $i$ 个方程的 $c_i={\displaystyle \frac{M}{m_i}}$
3. 计算 $c_i$ 模 $m_i$ 意义下的逆元 $c_i^{-1}$
4. 结果 $x=\sum _{i=1}^n{r}_i{c}_i{c}_i^{-1}(modM)$

由于一般题目中的 $m_1,m_2,...,m_r$ 不会交代是否两两互质，所以不能用中国剩余定理直接求解。

这是需要扩展中国剩余定理求解 。

扩展中国剩余定理(EXCRT)

问题改为求 $x$ 的最小非负整数解。

对于前两个方程 $x\equiv a_1(mod{m}_1)$ 和 $x\equiv a_2(mod{m}_2)$ ，

可以转化为丢番图方程 $x=m_{1}p+a_1=m_{2}q+a_2$ 。

所以，$m_{1}p-m_{2}q=a_2-a_1$。

由裴蜀定理

当 $(m_1,m_2)\nmid (a_2-a_1)$ 时，无解。

当 $(m_1,m_2)\mid (a_2-a_1)$ 时，有解。

由扩欧算法

$exgcd$ 函数返回值 $p,q$，可以得到特解 $p_0=p\times {\displaystyle \frac{a_2-a_1}{(m_1,m_2)}},q_0=q\times {\displaystyle \frac{a_2-a_1}{(m_1,m_2)}}$ ，

通解为 $P=p_0+{\displaystyle \frac{m_2}{(m_1,m_2)}}$，$Q=q_0-{\displaystyle \frac{m_1}{(m_1,m_2)}}$。

所以 $x=m_{1}P+a_1={\displaystyle \frac{m_1{m}_2}{(m_1,m_2)}}\times k+m_1\cdot p_0+a_1$。

数论中我们记 $a$ 和 $b$ 的最小公倍数为 $[a,b]$ 。

前两个方程等价合并为一个方程 $x\equiv a(modm)$，

其中 $a=m_1{p}_0+a_1$，$m=lcm(m_1,m_2)=[m_1,m_2]$。

所以 $n$ 个同余方程只要合并 $n-1$ 次，即可求解。

[AcWing204] 表达整数的奇怪方式

题目描述

给定 $2n$ 个整数 $a_1,a_2,\dots ,a_n$ 和 $m_1,m_2,\dots ,m_n$，求一个最小的非负整数 $x$，满足 $\mathrm{\forall }i\in [1,n],x\equiv m_i(mod{a}_i)$。

输入格式

第 $1$ 行包含整数 $n$。

第 $2\dots n+1$ 行：第 $i+1$ 行包含两个整数 $a_i$ 和 $m_i$，数之间用空格隔开。

输出格式

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

数据范围

$1\le a_i\le 2^{31}-1$,
$0\le m_i<a_i$
$1\le n\le 25$
所有 $m_i$ 的最小公倍数在 $64$ 位有符号整数范围内。

输入样例：

2
8 7
11 9

输出样例：

31

---

算法

扩展中国剩余定理模版题。

C++ 代码

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int N = 30;
ll a[N],m[N];
int n;
 
ll exgcd(ll a, ll b, ll &x, ll &y)
{
    if (b == 0){
        x = 1, y = 0;
        return a;
    }
    ll d,x1,y1;
    d = exgcd(b, a % b, x1, y1);
    x = y1, y = x1 - a / b * y1;
    return d;
}
 
ll excrt(ll a[], ll m[])
{
    ll m1, m2, a1, a2, p, q;
    m1 = m[1], a1 = a[1];
    for (int i = 2; i <= n; i++){
        m2 = m[i], a2 = a[i];
        ll d = exgcd(m1, m2, p, q);// d=(m1,m2)，同时计算m1p+m2q=1的特解(q=-y)
        if ((a2 - a1) % d) return -1;
        p = p * (a2 - a1) / d;// 这是m1p+m2q=a2-a1特解
        // 保证m2/d为正数
        p = (p % abs(m2 / d) + abs(m2 / d)) % abs(m2 / d);// 由于p可能是整数，需要加余取余操作变为最小正整数
        a1 = m1 * p + a1;
        m1 = m1 * m2 / d;
    }
    return (a1 % m1 + m1) % m1;// 同上
}
 
int main()
{
    cin >> n;
    // 这里与题目相反，一般是m是余数
    for (int i = 1; i <= n; i++) cin >> m[i] >> a[i];
    cout << excrt(a, m) << endl;
    
    return 0;
}

P4777 【模板】扩展中国剩余定理（EXCRT）

题目描述

给定 $n$ 组非负整数 $a_i,b_i$ ，求解关于 $x$ 的方程组的最小非负整数解。

$$\{\begin{array}{l}x\equiv b_1(\mathrm{mod}{a}_1)\\ x\equiv b_2(\mathrm{mod}{a}_2)\\ \dots \\ x\equiv b_n(\mathrm{mod}{a}_n)\end{array}$$

输入格式

输入第一行包含整数 $n$。

接下来 $n$ 行，每行两个非负整数 $a_i,b_i$。

输出格式

输出一行，为满足条件的最小非负整数 $x$。

输入输出样例 #1

输入 #1

3
11 6
25 9
33 17

输出 #1

809

说明/提示

对于 $100\mathrm{\%}$ 的数据，$1\le n\le {10}^5$，$1\le b_i,a_i\le {10}^{12}$，保证所有 $a_i$ 的最小公倍数不超过 ${10}^{18}$。

请注意程序运行过程中进行乘法运算时结果可能有溢出的风险。

数据保证有解。

算法

同样为扩展中国剩余定理的板子题。

但是注意不要越界，函数 $mul$ 和 $get$ 就在解决这个问题。

$mul$ 函数可以快速计算 $a\times b\mathrm{\%}m$，所以利用它可以计算 $p=p_0\times |{\displaystyle \frac{(a2-a1)}{(m1,m2)}}|\mathrm{\%}{\displaystyle \frac{m2}{(m1,m2)}}$。

代码

#include <bits/stdc++.h>
using namespace std;
#define IOS ios::sync_with_stdio(false)
typedef long long ll;
const int N = 1e5 + 10;
ll a[N], b[N];
int n;
 
ll mul(ll a, ll b, ll m)// 计算a * b % m
{
    auto res = 0ll;
    while (b){
        if (b & 1) res = (res + a) % m;
        a = (a + a) % m;
        b >>= 1;
    }
    return res;
}
 
ll exgcd(ll a, ll b, ll &x, ll &y)
{
    if (b == 0){
        x = 1, y = 0;
        return a;
    }
    ll d = exgcd(b, a % b, y, x);
    y -= a / b * x;
    return d;
}
 
inline ll get(ll a, ll b)// 保证最小整数
{
    return (a % b + b) % b;
}
 
ll excrt(ll a[], ll m[])
{
    ll m1, m2, a1, a2, p, q;
    m1 = m[1], a1 = a[1];
    auto ans = 0ll;
    for (int i = 2; i <= n; i++){
        m2 = m[i], a2 = a[i];
        ll d = exgcd(m1, m2, p, q);
        // 合并成ap + bq = c
        ll a = m1, b = m2, c = get(a2 - a1, m2);
        if (c % d) return -1;
        p = mul(p, c / d, b / d);
        ans = a1 + m1 * p;
        m1 = m2 / d * m1;
        ans = get(ans, m1);
        a1 = ans;
    }
    return ans;
}
 
 
int main()
{
    IOS;
    cin >> n;
    for (int i = 1; i <= n; i++) cin >> a[i] >> b[i];
    cout << excrt(b, a) << endl;
    return 0;
}