Số Học đồng Dư - Viblo

I. Phép đồng dư thức cơ bản

Đồng dư thức là phép toán lấy số dư của số này khi chia cho số khác, kí hiệu là %\%%. Ví dụ: 5%2=15 \% 2=15%2=1, khi đó có thể viết là 5≡15 \equiv 15≡1 (mod(mod(mod 2)2)2).

Phép đồng dư thức có tính chất phân phối đối với phép cộng, phép nhân và phép trừ, cụ thể như sau:

• (a+b)(a + b)(a+b) %\%% ccc =[(a= [(a=[(a %\%% c)+(bc) + (bc)+(b %\%% c)]c)]c)] %\%% ccc.
• (a−b)(a - b)(a−b) %\%% ccc =[(a= [(a=[(a %\%% c)−(bc) - (bc)−(b %\%% c)]c)]c)] %\%% ccc.
• (a×b)(a \times b)(a×b) %\%% ccc =[(a= [(a=[(a %\%% c)×(bc) \times (bc)×(b %\%% c)]c)]c)] %\%% ccc.

Riêng đối với phép chia, chúng ta không có tính chất phân phối, mà phải sử dụng một lí thuyết là Nghịch đảo modulo.

II. Nghịch đảo modulo của một số 1. Định nghĩa

Như chúng ta đều biết ở chương trình Toán phổ thông, nghịch đảo của một số nguyên aaa (kí hiệu a−1a^{-1}a−1) là số thỏa mãn: a.a−1=1a.a^{-1}=1a.a−1=1.

Đối với nghịch đảo modulo, ta cũng có khái niệm tương tự, nhưng là xét trên tập số dư khi chia cho MMM. Nghịch đảo modulo MMM của một số aaa (cũng kí hiệu a−1a^{-1}a−1) là số nguyên thỏa mãn: a.a−1≡1 (moda.a^{-1}\equiv1\ (moda.a−1≡1 (mod M)M)M) (Nói cách khác, a−1a^{-1}a−1 chính là 1a\frac{1}{a}a1​ %\%% M)M)M). Lấy ví dụ, nếu ta chọn M=109+7,a=2M={10}^9+7, a=2M=109+7,a=2 thì a−1=500000004a^{-1}=500000004a−1=500000004.

Không phải lúc nào cũng tồn tại a−1a^{-1}a−1. Chỉ khi GCD(a,M)=1GCD(a, M)=1GCD(a,M)=1 thì mới tồn tại a−1a^{-1}a−1 là nghịch đảo modulo MMM của aaa. Để tính nghịch đảo modulo của một số, ta có thể sử dụng hai giải thuật: Giải thuật Euclid mở rộng hoặc dựa trên định lý Fermat nhỏ (áp dụng giải thuật chia để trị tính ab a^b\ %\ cab ).

2. Giải thuật Euclid mở rộng

GCD(A,B)\text{GCD}(A,B)GCD(A,B) có một tính chất là luôn có thể biểu diễn dưới dạng phương trình:

Ax+By=GCD(A,B) (1)Ax+By=\text{GCD}(A,B) \ (1) Ax+By=GCD(A,B) (1)

Giải thuật Euclid mở rộng sử dụng để tìm một cặp số nguyên (x,y)(x,y)(x,y) thỏa mãn phương trình trên, và còn dùng để tính nghịch đảo modulo (mình sẽ nói tới ở phần sau). Cặp số (x,y)(x,y)(x,y) có thể có giá trị âm, hoặc bằng 000 đều được. Dưới đây tôi sẽ trình bày giải thuật tìm GCD(A,B)\text{GCD}(A,B)GCD(A,B) và một cặp (x,y)(x,y)(x,y) thỏa mãn phương trình.

long long d, x, y; // UCLN và cặp nghiệm (x, y). void extended_euclid(long long A, long long B) { if (B == 0) { x = 1; y = 0; d = A; } else { extended_euclid(B, A % B); long long temp = x; x = y; y = temp – (A / B) * y; } } int main() { cin >> A >> B; extended_euclid(A, B); cout << d << ' ' << x << ‘ ‘ << y; return 0; }

Cơ chế của giải thuật như sau: Ban đầu chương trình thực thi giống giải thuật Euclid, tới khi B=0,B=0,B=0, khi đó ta sẽ có AAA là ước chung lớn nhất của AAA và B,B,B, sau đó đặt x=1,y=0x=1,y=0x=1,y=0. Bởi vì B=0B=0B=0 và hiện tại GCD(A,B)=A\text{GCD}(A,B)=AGCD(A,B)=A nên phương trình (1)(1)(1) trở thành:

A.1+0.0=AA.1+0.0=A A.1+0.0=A

Sau đó chương trình gọi lại các lệnh dưới lời gọi đệ quy để tìm ra xxx và yyy. Chứng minh như sau:

• Sau khi gọi đệ quy, phương trình ở bước tiếp theo của giải thuật là:

B.x+(A%B).y=GCD(A,B) (2)B.x+(A \% B).y=\text{GCD}(A,B) \ (2) B.x+(A%B).y=GCD(A,B) (2)

• Thay A % B=A−⌊AB⌋.BA \ \% \ B=A-\left \lfloor{\frac{A}{B}} \right \rfloor .BA % B=A−⌊BA​⌋.B, phương trình (2)(2)(2) trở thành:

B.x+(A-\left \lfloor{\frac{A}{B}} \right \rfloor.B).y=\text{GCD}(A,B)$$ $$\Leftrightarrow A.y+B.(x-\left \lfloor{\frac{A}{B}} \right \rfloor.y)=\text{GCD}(A,B)

• Từ đây được công thức đệ quy:

x′=y;y′=x−⌊AB⌋.yx' = y; y' = x - \left \lfloor{\frac{A}{B}} \right \rfloor.y x′=y;y′=x−⌊BA​⌋.y

• Như vậy từ bước cơ bản x=1,y=0;x=1,y=0;x=1,y=0; chương trình sẽ tiếp tục tính ngược lên để ra được x,yx,yx,y thỏa mãn phương trình ban đầu. Độ phức tạp giải thuật là O(log⁡(max⁡(A,B)))O\Big(\log\big(\text{max}⁡(A,B)\big)\Big)O(log(max⁡(A,B))).

Ngoài ra, giải thuật Euclid mở rộng còn có thể sử dụng để giải phương trình Diophantine tuyến tính, sẽ được đề cập tới ở một bài viết khác.

3. Tính toán nghịch đảo Modulo của một số

Sử dụng giải thuật Euclid mở rộng

Như đã trình bày ở trên, theo giải thuật Euclid mở rộng, nếu GCD(a,M)=1GCD\left(a,M\right)=1GCD(a,M)=1, ta luôn tìm được xxx và yyy thỏa mãn: a.x+M.y=1a.x+M.y=1a.x+M.y=1. Mà M.yM.yM.y chia hết cho MMM, do đó phương trình trở thành:

a.x≡1(mod M)a.x \equiv 1 (\text{mod} \ M) a.x≡1(mod M)

Từ đây suy ra xxx chính là a−1a^{-1}a−1. Tuy nhiên trong giải thuật Euclid mở rộng, xxx có thể mang giá trị âm, nên ta sẽ điều chỉnh một chút để tính được giá trị a−1a^{-1}a−1 luôn không âm. Code dưới đây sẽ tái sử dụng lại đoạn code giải thuật Euclid mở rộng ở phía trên:

long long modulo_inverse(long long a, long long M) { extended_euclid(a, M); // a và M không nguyên tố cùng nhau, không tồn tại nghịch đảo modulo M của a. if (d != 1) return -1; // Do x có thể âm, ta làm dương nó. return (x % M + M) % M; }

Tính nghịch đảo modulo bằng định lý Fermat nhỏ

Theo định lý Fermat nhỏ, ta có: Nếu MMM là số nguyên tố và aaa không chia hết cho MMM thì:

a^{M-1}\equiv 1 \ (\text{mod} \ M)$$ hay nói cách khác: $$a\times a^{M-2} \equiv 1 \ (\text{mod} \ M)

Điều này tương đương với việc nếu MMM là số nguyên tố và aaa không chia hết cho MMM thì aM−2a^{M-2}aM−2 chính là nghịch đảo modulo MMM của aaa, cũng tương đương với aM−2a^{M-2}aM−2 %\%% MMM là nghịch đảo modulo MMM của aaa.

long long power_mod(long long a, long long b, long long M) // Tính a^b % M. { if (b == 0) return 1; if (b == 1) return a; long long half = power_mod(a, b / 2, M) % M; if (b % 2 == 0) return (half * half) % M; else return (((half * half) % M) * a) % M; } long long modulo_inverse(int a, int M) { return power_mod(a, M – 2, M); }

4. Áp dụng nghịch đảo modulo để tính ab % c\frac{a}{b} \ \% \ cba​ % c

Mình đã đề cập ở mục 111, phép chia không có tính chất phân phối đối với phép đồng dư thức giống như các phép cộng, trừ và nhân. Để tính ab % c,\frac{a}{b} \ \% \ c,ba​ % c, ta làm như sau:

• Tách ab=(a×1b) % c=(a×b−1) % c,\frac{a}{b} = \left(a \times \frac{1}{b}\right) \ \% \ c =\left(a \times b^{-1}\right) \ \% \ c,ba​=(a×b1​) % c=(a×b−1) % c, trong đó b−1b^{-1}b−1 là nghịch đảo modulo ccc của bbb.
• Sau đó áp dụng tính chất phân phối của phép nhân đối với phép đồng dư thức, lúc này phép chia modulo trở thành phép nhân với nghịch đảo modulo. Lưu ý, tùy vào giá trị ccc mà ta chọn cách tìm nghịch đảo modulo thích hợp (ccc có là số nguyên tố hay không).

Cài đặt:

long long modulo_divide(long long a, long long b, long long c) { long long inverse = modulo_inverse(b, c); return (a % c * inverse) % c; } III. Một số kiến thức nâng cao khác

1. Bậc lũy thừa theo modulo NNN (Multiplicative Order)

Xét hai số nguyên aaa và NNN nguyên tố cùng nhau, bậc lũy thừa của aaa theo modulo NNN là số nguyên dương KKK nhỏ nhất thỏa mãn: aK≡1 (mod N)a^K \equiv 1 \text{ } (mod \text{ } N)aK≡1 (mod N), kí hiệu là ordN(a)ord_N(a)ordN​(a).

Theo định lý Euler, vì aaa và NNN là hai số nguyên tố cùng nhau nên aϕ(N)≡1 (mod N),a^{\phi(N)} \equiv 1 \ (\text{mod} \ N),aϕ(N)≡1 (mod N), với ϕ(N)\phi(N)ϕ(N) là số lượng số nguyên dương không vượt quá NNN và nguyên tố cùng nhau với NNN. Mà ϕ(N)≤N\phi(N) \le Nϕ(N)≤N, do đó ordN(a)≤Nord_N(a) \le NordN​(a)≤N, vậy để tìm ordN(a)ord_N(a)ordN​(a) chỉ cần duyệt một vòng lặp từ 111 tới NNN với độ phức tạp O(N−1)O(N - 1)O(N−1).

int find_m_order(int a, int N) { int mul = 1; for (int i = 1; i <= N; ++i) { mul = (mul * a) % N; if (mul == 1) return i; } }

2. Tiêu chuẩn Euler (Euler's Criterion)

Đầu tiên, ta làm quen với khái niệm Thặng dư bậc hai: Một số nguyên qqq được gọi là thặng dư bậc hai theo modulo NNN nếu như nó đồng dư với một số chính phương theo modulo N,N,N, nghĩa là tồn tại một số nguyên xxx sao cho x2≡q (mod N)x^2 \equiv q \ (\text{mod} \ N)x2≡q (mod N).

Trong lý thuyết số, tiêu chuẩn Euler là một công thức dùng để xác định xem một số nguyên có phải là một thặng dư bậc hai theo modulo PPP (với PPP là một số nguyên tố) hay không. Theo đó, xét hai số nguyên aaa và PPP nguyên tố cùng nhau, trong đó PPP là một số nguyên tố lẻ. Ta có công thức sau:

aP−12≡{1 (mod P),neˆˊu a laˋ thặng dư bậc hai của P.−1 (mod P),neˆˊu a khoˆng laˋ thặng dư bậc hai của P.a^{\frac{P - 1}{2}} \equiv \begin{cases}1 \text{ }(\text{mod} \ P),&\text{nếu }a \text{ là thặng dư bậc hai của }P. \\ -1\text{ }(\text{mod} \ P),& \text{nếu }a \text{ không là thặng dư bậc hai của }P.\end{cases} a2P−1​≡{1 (mod P),−1 (mod P),​neˆˊu a laˋ thặng dư bậc hai của P.neˆˊu a khoˆng laˋ thặng dư bậc hai của P.​

Đối với trường hợp P=2,P=2,P=2, mọi số nguyên đều là thặng dư bậc hai theo modulo PPP.

Ví dụ, xét P=7P = 7P=7, ta có a=2a = 2a=2 là thặng dư bậc hai của 777, vì tồn tại hai số nguyên x=3x = 3x=3 và x=4x = 4x=4 thỏa mãn a≡x2 (mod P)a \equiv x^2 \text{ } (mod \text{ } P)a≡x2 (mod P).

long long power_mod(long long a, long long b, long long P) { if (b == 0) return 1; if (b == 1) return a; long long half = power_mod(a, b / 2, P) % P; if (b % 2 == 0) return (half * half) % P; else return (((half * half) % P) * (a % P)) % P; } // Kiểm tra N có phải thặng dư bậc 2 của P hay không. bool check_quadratic_residue(long long N, long long P) { if (P == 2) return true; else return (power_mod(N, (P - 1) / 2, P) == 1); }

Trong trường hợp NNN và PPP cùng có dạng 4i+3 (i>0)4i + 3 \ (i > 0)4i+3 (i>0), thì giá trị xxx thỏa mãn x2≡N (mod P)x^2 \equiv N \ (\text{mod} \ P)x2≡N (mod P) (nếu tồn tại) chỉ có thể là: x=± NP+14x=\pm \text{ }N^{\frac{P + 1}{4}}x=± N4P+1​. Dựa vào nhận xét này ta có thể tính nhanh ra giá trị xxx. Chứng minh nhận xét như sau:

• Theo định lý Euler, ta có: NP−12 % P=1N^{\frac{P - 1}{2}} \ \% \ P = 1N2P−1​ % P=1.
• Nhân cả hai vế với NNN:

NP+12 % P=N % P (1)N^{\frac{P + 1}{2}} \ \% \ P = N \ \% \ P \ (1) N2P+1​ % P=N % P (1)

• Lại có: x2≡N (mod P)x^2 \equiv N \ (\text{mod} \ P)x2≡N (mod P).       ⇔x2≡NP+12 (mod P)\Leftrightarrow x^2 \equiv N^{\frac{P+1}{2}} \ (\text{mod} \ P)⇔x2≡N2P+1​ (mod P) (do đẳng thức (1)(1)(1)).       ⇔x2≡N2i+2 (mod P)\Leftrightarrow x^2 \equiv N^{2i + 2} \ (\text{mod} \ P)⇔x2≡N2i+2 (mod P) (do N=4i+3N=4i + 3N=4i+3).       ⇔x≡ Ni+1 (mod P)\Leftrightarrow x \equiv \ N^{i + 1} \ (\text{mod} \ P)⇔x≡ Ni+1 (mod P).       ⇔x≡± NP+14 (mod P)\Leftrightarrow x \equiv \pm \ N^{\frac{P + 1}{4}} \ (\text{mod} \ P)⇔x≡± N4P+1​ (mod P) (do P=4i+3P=4i + 3P=4i+3).

Cài đặt:

int find_quadratic_residue(int N, int P) { // P và N không ở đúng dạng 4i + 3, không tính được theo cách này. if (P % 4 != 3 || N % 4 != 3) return -1; int x = power_mod(N, (P + 1) / 2, P); // Kiểm tra giá trị x thứ nhất có hợp lệ không. if ((x * x) % P == N % P) return x; // Kiểm tra giá trị x thứ hai có hợp lệ không. x = P - x; if ((x * x) % P == N % P) return x; // Nếu không tồn tại x, kết luận N không phải thặng dư bậc 2 của P return -1; } IV. Tài liệu tham khảo

• https://en.wikipedia.org/wiki/Quadratic_residue
• https://en.wikipedia.org/wiki/Euler%27s_criterion
• https://vnoi.info/wiki/translate/he/So-hoc-Phan-1-Modulo-gcd.md
• https://vnoi.info/wiki/algo/math/modular-inverse.md
• https://en.wikipedia.org/wiki/Modular_multiplicative_inverse
• https://en.wikipedia.org/wiki/Multiplicative_order
• https://vi.wikipedia.org/wiki/%C4%90%E1%BB%8Bnh_l%C3%BD_Euler
• https://www.geeksforgeeks.org/find-square-root-under-modulo-p-set-1-when-p-is-in-form-of-4i-3/

©️ Tác giả: Vũ Quế Lâm từ Viblo