\begin{cases} x\equiv a_1\ \left( \text{mod\ }p_1 \right) \\ x\equiv a_2\ \left( \text{mod\ }p_2 \right) \\ x\equiv a_3\ \left( \text{mod\ }p_3 \right) \end{cases}
以前做 CRT 一直是用那个构造性的暴力式子 ( q_{ji} 表示 p_j 在 p_i 下的逆元 ):
x \equiv \sum_{i=1}^{n}{( a_i \cdot \prod_{j \ne i}{(p_j \cdot q_{ji})} )} \mod {(\prod_{i=1}^{n}{p_i})}
(不过这个式子无比好记&&好用)
现在我们来考虑用 exgcd 去合并 CRT 方程。
直接先考虑俩方程
x \equiv a_1 \mod p_1
x \equiv a_2 \mod p_2
搞成某种加法形式
x + y_1p_1 = a_1
x + y_2p_2 = a_2
减一下
y_1p_1-y_2p_2=a_1-a_2
我习惯把减号改成加号
y_1p_1+y_2p_2=a_1-a_2
直接暴力 exgcd 这个方程,无解的条件是当且仅当 \gcd(p_1,p_2) \nmid (a_1-a_2) 。
求出解后记得 \mod p_2 ,这相当于合并方程了(注意 CRT 方程组最后的解的模数是 \prod_{i=1}^{n}{p_i} 。)
然后回代出第一个方程的特解 x_1
x_1 = a_1-y_1p_1
(代码上的一个小技巧:因为这特解和 a_1 是相等同余的并且我们还要合并其他方程,所以直接更新 a_1 )
然后我们考虑通解是这玩意 :
x \equiv x_1 \mod \mathrm{lcm}(a_1,a_2)
所以记得变一下模数。(可以维持那个特解不用变,因为在 \mathrm{lcm} 下是同余的)
大概的(伪)代码:
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 |
int P=p1,A=a1; for (i=2 to n) { exgcd(P,pi,x,y); G=gcd(P,pi); x= ((A-ai)/G)*x %pi; // exgcd for Diophantus eqution, mind that we need to modulo pi // No sol if ((P-pi)%G!=0) return & No sol; A-=x*P; // solve x_1, mind that x_1 equiv to A P=P/G*pi; // lcm is the new modulo of new eqution A%=P; // in case of overflow } return (A%P+P)%P; |