【学习笔记】CRT & ExtendedCRT

中国剩余定理是用于求解线性同余方程组模数两两互质的情形的一种方法。

说人话就是求解以下线性同余方程组：

$$\begin{cases} x \equiv a_1 & \pmod {m_1} \\ x \equiv a_2 & \pmod {m_2} \\ & \dots \\ x \equiv a_n & \pmod {m_n} \end{cases}$$

其中 $m_1, m_2\cdots m_n$ 两两互质。

其实中国剩余定理只要背结论就好了

　$\text{CRT}$

我们令 $m = \prod_{i=1}^nm_i$ , $M_i = \frac{m}{m_i}$ , $t_i$ 为同余方程 $M_it_i \equiv 1 \pmod{m_i}$ 的一个解

那么上述线性同余方程组的一个解为

$$x = \sum_{i=1}^n a_iM_it_i$$

证明什么的就略过了

$\text{ExtendedCRT}$

$\text{ExtendedCRT}$ 主要用于处理模数不互质的情况，处理方式有两种，一种是拆一个同余方程为多个同余方程，其中每个的模数两两互质，另一种方法是考虑合并同余方程组。这里讲第二种方法（其实是更普遍的方法）。

这其实就是一个子问题，只要我们能处理 $n = 2$ 的情形，这就意味着我们合并了两个同余方程，那么后面的结果就可以归纳得出。

考虑 $n=2$ 的情形

$$\begin{cases} x \equiv a_1 & \pmod {m_1} \\ x \equiv a_2 & \pmod {m_2} \end{cases}$$

我们有

$$\begin{aligned} x = a_1+\lambda m_1 &= a_2+\mu m_2 \\ \lambda m_1-\mu m_2 &= a_2-a_1 \end{aligned}$$

我们令 $m_0 = (m_1, m_2)$

• 若 $a_1 \not\equiv a_2 \pmod{m_0}$ ，则方程无解
• 若 $a_1 \equiv a_2 \pmod {m_0}$ ，则原方程可化为

$$\frac{m_1}{m_0}\lambda-\frac{m_2}{m_0}\mu=\frac{a_2-a_1}{m_0}$$

先令右式等于 $1$ ，求得特解，再将解乘上 $\frac{a_2-a_1}{m_0}$ ，得到两个解 $\lambda_0$, $\mu_0$ , 求得通解

$$\begin{aligned} \lambda &= \lambda_0 - \frac{m_2}{m_0}t \\ \mu &= \mu_0 + \frac{m_1}{m_0}t \\ t & \in \Z \end{aligned}$$

带入关于 $x$ 的方程， 得

$$\begin{aligned} x &= a_1+\lambda m_1 \\ &= a_1 + (\lambda_0-\frac{m_2}{m_0}t)m_1 \\ &= a_1 + \lambda_0m_1 - \frac{m_1m_2}{m_0}t \\ &= a_1 + \lambda_0m_1 - [m_1,m_2]t \\ x & \equiv a_1+\lambda_0m_1 \pmod{[m_1, m_2]} \end{aligned}$$

这样我们就完成了两个同余方程的合并，按照这种步骤全部合并就会得到解了。

$\text{Code}$

// ------------------------ CRT -----------------
inline int CRT() {
    int ans = 0, lcm = 1;
    for (int i = 1; i <= n; i++) lcm *= m[i];
    for (int i = 1; i <= n; i++) {
        int tmp = lcm / m[i];
        (ans += a[i] * tmp * ModReverse(tmp, m[i])) %= lcm;
    }
    return (ans + lcm) % lcm;
}

// ----------------- Extended CRT --------------

inline int ExtendedCRT() {
	for (int i = 2, x, y, z; i <= n; i++) {
		z = ExtendedGCD(m[1], m[i], x, y);
		if ((a[1] - a[i]) % z != 0) return -1;
		ExtendedGCD(m[1] / z, m[i] / z, x, y);
		x *= (a[i] - a[1]) / z, y *= (a[1] - a[i]) / z;
		a[1] += x * m[1], (m[1] /= z) *= m[i], a[1] %= m[1];
		if (a[1] < 0) a[1] += m[1];
	}
	return a[1];	
}