1. 程式人生 > 其它 >中國剩餘定理簡析(python實現)

中國剩餘定理簡析(python實現)

簡述中國剩餘定理CRT,使用python進行演算法實現

中國剩餘定理CRT

正整數m1,m2,...,mk兩兩互素,對b1,b2,...,bk的同餘式組為

\[\begin{cases} x \equiv b_1\; mod \;m_1\\ x \equiv b_2\; mod \;m_2\\ \quad\quad\vdots\\ x \equiv b_k\; mod \;m_k\\ \end{cases} \]

在mod M

\[M = \prod_{i = 1}^{k}m_i \]

的情況下有唯一解

\[x = (\sum_{i=1}^k b_iM_iM'_i)\;mod\;M \]

其中

\[M_i = \frac{M}{m_i} \]\[M'_i = M_i^{-1}\;mod\;m_i \]

python程式碼實現:

import gmpy2

def crt(b,m):
    #判斷是否互素
    for i in range(len(m)):
        for j in range(i+1,len(m)):
            if gmpy2.gcd(m[i],m[j]) != 1:
                print("m中含有不是互餘的數")
                return -1
    #乘積
    M = 1
    for i in range(len(m)):
        M *= m[i]
    #求M/mi
    Mm = []
    for i in range(len(m)):
        Mm.append(M // m[i])
    #求Mm[i]的乘法逆元
    Mm_ = []
    for i in range(len(m)):
        _,a,_ = gmpy2.gcdext(Mm[i],m[i])
        Mm_.append(int(a % m[i]))
    #求MiM'ibi的累加
    y = 0
    for i in range(len(m)):
        print(Mm[i] * Mm_[i] * b[i])
        y += (Mm[i] * Mm_[i] * b[i])
    y = y % M
    return y