I. Giới thiệu chung

Nếu như là một học sinh chuyên Toán, chắc hẳn các bạn đã nghe đến khái niệm Phương trình Diophantine tuyến tính (Linear Diophantine Equation). Đó là phương trình bậc nhất hai ẩn có dạng:

$$ax + by = c$$

với $a, b, c$ là các số nguyên cho trước.

Giải phương trình Diophantine tuyến tính là một vấn đề không đơn giản, nó yêu cầu người giải phải có kiến thức khá tốt về Toán học, cũng như một số kiến thức bổ trợ mà thường được nói đến nhiều hơn trong Tin học như giải thuật Euclid và giải thuật Euclid mở rộng. Trong chuyên đề này, chúng ta sẽ cùng nhau nghiên cứu cách giải phương trình trên và một số bài tập dạng biến thể của nó như:

• Tìm một nghiệm của phương trình.
• Tìm tất cả các nghiệm của phương trình.
• Đếm số nghiệm của phương trình và số nghiệm của phương trình trên một đoạn ràng buộc.
• Tìm nghiệm của phương trình thỏa mãn $x + y$ đạt giá trị nhỏ nhất.

Trước khi nghiên cứu bài viết này, các bạn cần có kiến thức về Giải thuật Euclid mở rộng. Nếu như chưa nắm được kiến thức này, các bạn có thể nghiên cứu lại về nó tại đây.

II. Giải thuật cơ sở tìm một nghiệm của phương trình

Phân tích toán học

Trước hết, chúng ta sẽ cùng nghiên cứu cách tìm ra một nghiệm của phương trình. Để cho đơn giản, tôi sẽ gọi nghiệm là nghiệm riêng $(x_0, y_0)$ của phương trình. Từ nghiệm riêng này, chúng ta có thể tìm ra mọi nghiệm của phương trình Diophantine tuyến tính.

Trường hợp suy biến của phương trình là $a = 0, b = 0$. Dễ dàng nhận ra phương trình sẽ có vô số nghiệm hoặc vô nghiệm, tùy thuộc vào $c = 0$ hay $c \ne 0$. Trường hợp này sẽ tạm thời được bỏ qua trong bài viết, vì nó quá đơn giản.

Khi $a \ne 0, b \ne 0;$ phương trình $ax + by = c$ có thể được biến đổi thành một trong hai phương trình dưới đây đều được:

$$\begin{cases}ax \equiv c \ (\text{mod } b) \ (1) \\ by \equiv c \ (\text{mod } a) \ (2) \end{cases}$$

Không mất tính tổng quát, hãy giả sử $b \ne 0$ và biến đổi phương trình ban đầu về phương trình $(1)$. Khi $a$ và $b$ là hai số nguyên tố cùng nhau, thì nghiệm của phương trình sẽ có dạng:

$$\begin{cases}x_0 \equiv c.a^{-1} \ (\text{mod }b) \\ y_0 = \frac{c - ax_0}{b} \end{cases} (3)$$

với $a^{-1}$ là nghịch đảo modulo $b$ của $a$

Còn nếu như $a$ và $b$ không nguyên tố cùng nhau, thì gọi $g = \text{GCD}(a, b)$. Khi đó $a$ và $b$ đều chia hết cho $g,$ vì thế phương trình sẽ chỉ có nghiệm nếu như $c$ cũng chia hết cho $g$. Nghiệm $x_0$ của nó sẽ có dạng:

$$\frac{a}{g} \times x_0 \equiv \frac{c}{g} \ \left(\text{mod } \frac{b}{g}\right) \ (4)$$

Bởi vì $g = \text{GCD}(a, b);$ nên $\frac{a}{g}$ và $\frac{b}{g}$ sẽ lại là hai số nguyên tố cùng nhau. Phương trình $(4)$ sẽ quay trở lại dạng $(3)$ và nghiệm của nó sẽ là:

$$\begin{cases}x_0 \equiv \frac{c}{g} \times \left(\frac{a}{g}\right)^{-1} \ \left(\text{mod } \frac{b}{g}\right). \\ y_0 = \frac{c - ax_0}{b} \cdot \end{cases}$$

Như đã nói, ta sẽ gọi đây là nghiệm riêng của phương trình.

Ý tưởng giải thuật

Phân tích toán học là vậy, còn khi giải quyết bằng lập trình, chúng ta sẽ sử dụng giải thuật Euclid mở rộng. Trước hết, giả sử $a$ và $b$ đều là các số nguyên không âm. Khi sử dụng giải thuật Euclid mở rộng với hai số này, ta sẽ thu được ước chung lớn nhất $g$ của chúng và hai số $x_g, y_g$ thỏa mãn:

$$ax_g + by_g = g \ (5)$$

Nếu $c$ chia hết cho $g,$ thì phương trình Diophantine tuyến tính sẽ có nghiệm (như đã chứng minh ở trên). Nhân cả hai vế của phương trình $(5)$ với $\frac{c}{g},$ ta có:

$$a \times x_g \times \frac{c}{g} + b \times y_g \times \frac{c}{g} = c$$

Vì thế, nghiệm riêng của phương trình Diophantine sẽ là:

$$\begin{cases}x_0 = x_g \times \frac{c}{g} \\ y_0 = y_g \times \frac{c}{g} \end{cases}$$

Giải thuật trên vẫn hoàn toàn chính xác khi $a, b$ hoặc cả hai số cùng là các số nguyên âm. Ta chỉ cần thay đổi dấu của $x_0$ hoặc $y_0$ sau khi tìm ra nghiệm nếu như $a < 0$ hoặc $b < 0$ tương ứng.

Cài đặt

// Trả về gcd(a, b) và cập nhật nghiệm (x_g, y_g) của phương trình.
int extended_euclid(int a, int b, int &x, int &y)
{
    if (b == 0)
    {
        x = 1;
        y = 0;
        return a;
    }
	
    int x1, y1;
    int d = extended_euclid(b, a % b, x1, y1);
    x = y1;
    y = x1 - y1 * (a / b);
	
    return d;
}

// Trả về true và cập nhật nghiệm (x0, y0) nếu có nghiệm.
// Ngược lại trả về false và không cập nhật nghiệm.
bool find_one_solution(int a, int b, int c, int &x0, int &y0, int g)
{
    g = extended_euclid(abs(a), abs(b), x0, y0);
    
    if (c % g != 0)
        return false;
		
    x0 *= c / g;
    y0 *= c / g;
	
    if (a < 0)
        x0 = -x0;
    if (b < 0)
        y0 = -y0;
		
    return true;
}

void solve_diophantine_equation(int a, int b, int c)
{
    // Trường hợp suy biến: a = 0 và b = 0.
    if (a == 0 && b == 0)
    {
        if (c == 0)
            cout << "The equation has infinite solutions";
        else 
            cout << "Solution is not existed";
    }
	
    // Giải phương trình và tìm một nghiệm hoặc thông báo vô nghiệm.
    int x0, y0, g;
    bool has_solution = find_one_solution(a, b, c, x0, y0, g);
	
    if (has_solution)
        cout << x0 << ' ' << y0;
    else 
        cout << "Solution is not existed";
}

Đánh giá

Độ phức tạp của giải thuật bằng với độ phức tạp của thuật toán Euclid mở rộng, do đó độ phức tạp tổng quát là: $O\Big(\log\big(\text{max}(a, b)\big)\Big)$.

III. Một số vấn đề mở rộng

1. Tìm tất cả các nghiệm của phương trình Diophantine

Ý tưởng

Từ việc tìm ra nghiệm riêng $(x_0, y_0)$ của phương trình, ta có thể tìm ra tất cả các nghiệm của phương trình.

Đặt $g = \text{GCD}(a, b)$ và giả sử đã tìm ra nghiệm riêng $(x_0, y_0)$ thỏa mãn phương trình:

$$a \times x_0 + b \times y_0 = c$$

Ta sẽ biến đổi phương trình bằng cách cộng thêm $\frac{b}{g}$ vào $x_0$ và trừ đi $\frac{a}{g}$ ở $y_0;$ khi đó tính cân bằng của phương trình vẫn không bị phá vỡ:

$$a \times \left(x_0 + \frac{b}{g}\right) + b \times \left(y_0 - \frac{a}{g}\right) = a \times x_0 + b \times y_0 + a \times \frac{b}{g} - b \times \frac{a}{g} = c$$

Cách biến đổi trên có thể lặp đi lặp lại vô hạn lần, mỗi lần đều cho ra một phương trình mới với tính cân bằng vẫn giữ nguyên. Vì thế, ta sẽ có các nghiệm của phương trình ở dạng tổng quát là:

$$\begin{cases}x = x_0 + k \times \frac{b}{g} \\ y = y_0 - k \times \frac{a}{g} \end{cases}$$

Tất cả các cặp $(x, y)$ có dạng này sẽ chính là tập toàn bộ các nghiệm của phương trình Diophantine tuyến tính hai ẩn.

2. Tìm mọi nghiệm của phương trình trong một đoạn

Ý tưởng

Dễ nhận thấy, nếu như không có ràng buộc gì thêm về nghiệm, thì từ những gì đã chứng minh ở phần trước, ta sẽ có vô số nghiệm của phương trình Diophantine. Do đó, thông thường chúng ta sẽ tìm nghiệm của phương trình trên một đoạn nào đó.

Giả sử ta có hai đoạn $[xmin, xmax]$ và $[ymin, ymax]$ là ràng buộc của nghiệm $x$ và $y$ chúng ta muốn tìm.

Trường hợp dễ nhận thấy là nếu $a = 0$ hoặc $b = 0$ thì phương trình chỉ có duy nhất một nghiệm và ta sẽ biết ngay nó có thuộc đoạn ràng buộc hay không. Do đó ta không xét trường hợp này ở đây.

Trước tiên, tìm ra nghiệm riêng $(x_0, y_0)$ bất kì của phương trình Diophantine. Sau đó, ta sẽ tìm hai giá trị $x_1, x_2$ từ $x_0$ theo công thức nghiệm tổng quát sao cho:

$$\begin{cases}x_1 \ge xmin; & x_1 \text{ nhỏ nhất} \\ x_2 \le xmax; & x_2 \text{ lớn nhất} \end{cases}$$

Đặt tập hợp $A = \{x \ | \ x \in [x_1, x_2]\}$.

Hoàn toàn tương tự, ta có thể tìm ra hai giá trị $(y_1, y_2)$ từ $y_0$ theo công thức nghiệm tổng quát sao cho:

$$\begin{cases}y_1 \ge ymin; & y_1 \text{ nhỏ nhất} \\ y_2 \le ymax; & y_2 \text{ lớn nhất} \end{cases}$$

Kế đến, tìm giá trị $x$ tương ứng với $y_1$ và $y_2$ vừa tìm ra theo công thức nghiệm $y = \frac{c - ax}{b}$. Ta gọi hai giá trị $x$ tương ứng trong các nghiệm này là $x_1'$ và $x_2'$. Gọi tập hợp $B = \{x \ | \ x \in [x_1', x_2']\}$.

Số lượng nghiệm của phương trình thỏa mãn $x \in [xmin, xmax]$ và $y \in [ymin, ymax]$ khi đó sẽ là số lượng phần tử của tập giao giữa hai tập hợp $A$ và $B$:

$$\text{D}(xmin, xmax, ymin, ymax) = |A \cap B|$$

Cài đặt

Cài đặt dưới đây sẽ đếm số nghiệm $(x, y)$ của phương trình Diophantine thỏa mãn:

$$\begin{cases}x \in [xmin, xmax] \\ y \in [ymin, ymax] \end{cases}$$

Lưu ý rằng, trước hết ta sẽ biến đổi phương trình về dạng $a, b$ là hai số nguyên tố cùng nhau để đơn giản hóa. Điều này được thực hiện bằng cách chia cả $a$ và $b$ cho ước chung lớn nhất $g$ của nó, bởi vì ta đều biết:

$$ax + by = c \Leftrightarrow \frac{a}{g}x + \frac{b}{g}y = \frac{c}{g}$$

Sau khi đã chia $a$ và $b$ cho $g,$ thì công thức tổng quát của nghiệm sẽ chỉ còn là:

$$\begin{cases} x = x_0 + k \times b \\ y = y_0 - k \times a\end{cases}$$

Cài đặt cũng sẽ sử dụng lại hàm tìm nghiệm riêng của phương trình find_one_solution().

// Sử dụng công thức tổng quát để tìm ra nghiệm mới từ nghiệm riêng.
void shift_solution(int &x, int &y, int a, int b, int k)
{
    x += k * b;
    y -= k * a;
}

int all_solutions_in_range(int a, int b, int c, 
                           int xmin, int xmax, int ymin, int ymax)
{
    int x, y, g;

    // Nếu không có nghiệm riêng thì cũng không có nghiệm trong đoạn ràng buộc.
    if (!find_one_solution(a, b, c, x, y, g))
        return 0;

    // Đưa a và b thành hai số nguyên tố cùng nhau.
    a /= g;
    b /= g;
	
    // Lưu dấu của a và b vào hai biến sign_a và sign_b.
    int sign_a = a > 0 ? 1 : -1;
    int sign_b = b > 0 ? 1 : -1;
	
    // Tìm nghiệm lx1 nhỏ nhất thỏa mãn xmin <= lx1 <= xmax. 
    shift_solution(x, y, a, b, (xmin - x) / b);
    if (x < xmin)
        shift_solution(x, y, a, b, sign_b);
    if (x > xmax)
        return 0;
    int lx1 = x;
	
    // Tìm nghiệm rx1 lớn nhất thỏa mãn xmin <= rx1 <= xmax.
    shift_solution(x, y, a, b, (xmax - x) / b);
    if (x > xmax)
        shift_solution(x, y, a, b, -sign_b);
    int rx1 = x;

    // Tìm nghiệm y nhỏ nhất thuộc đoạn [ymin, ymax]. 
    // Sau đó tìm nghiệm lx2 tương ứng với nó.
    shift_solution(x, y, a, b, -(miny - y) / a);
    if (y < ymin)
        shift_solution(x, y, a, b, -sign_a);
    if (y > ymax)
        return 0;
    int lx2 = x;

    // Tìm nghiệm y lớn nhất thuộc đoạn [ymin, ymax]. 
    // Sau đó tìm nghiệm rx2 tương ứng với nó.
    shift_solution(x, y, a, b, -(ymax - y) / a);
    if (y > ymax)
        shift_solution(x, y, a, b, sign_a);
    int rx2 = x;
	
    // Tìm đoạn [lx, rx] là giao của hai đoạn [lx1, rx1] và [lx2, rx2].
    if (lx2 > rx2)
        swap(lx2, rx2);
    int lx = max(lx1, lx2);
    int rx = min(rx1, rx2);
	
    // Kết luận số nghiệm của phương trình là độ dài đoạn [lx, rx].
    if (lx > rx)
        return 0;
    return (rx - lx) / abs(b) + 1;
}

Một khi đã tìm ra đoạn nghiệm $[l_x, r_x]$ cua giá trị $x$ thì ta cũng dễ dàng tìm được cụ thể các nghiệm đó theo công thức:

$$\begin{cases}x = l_x + k \times \frac{b}{g}; \forall k \ge 0 \text{ và } x \le r_x \\ y = \frac{c - ax}{b}\end{cases}$$

3. Tìm nghiệm sao cho $x + y$ đạt giá trị nhỏ nhất

Đôi khi, người ra đề có thể làm khó chúng ta bằng cách yêu cầu tìm nghiệm $(x, y)$ thỏa mãn $(x + y)$ đạt giá trị nhỏ nhất và $x \in [xmin, xmax]; y \in [ymin, ymax]$. Khi đó, ta sẽ thực hiện giống như ý tưởng tìm nghiệm thuộc một đoạn ràng buộc: Điều chỉnh nghiệm từ nghiệm riêng để thu được các nghiệm thỏa mãn điều kiện mong muốn.

Theo công thức nghiệm tổng quát của phương trình, ta biết rằng:

$$\begin{cases}x = x_0 + k \times \frac{b}{g} \\ y = y_0 - k \times \frac{a}{g}\end{cases}$$

Khi đó, $x + y$ sẽ trở thành:

$$x + y = x_0 + y_0 + k \times \left(\frac{b}{g} - \frac{a}{g}\right) = x + y + k \times \frac{b - a}{g}$$

Do đó, nếu muốn $x + y$ đạt giá trị nhỏ nhất, thì ta cần xử lý ba trường hợp:

• Nếu $a > b,$ cần tìm giá trị $k$ nhỏ nhất sao cho $x, y$ vẫn thuộc đoạn ràng buộc.
• Nếu $a < b,$ cần tìm giá trị $k$ lớn nhất sao cho $x, y$ vẫn thuộc đoạn ràng buộc.
• Nếu $a = b$ thì mọi cặp nghiệm $(x, y)$ của phương trình đều có tổng $x + y$ bằng nhau.

Cài đặt cụ thể xin dành lại cho bạn đọc. Các bạn có thể tham khảo một số bài tập vận dụng phương trình Diophantine ở dưới đây:

• Spoj - Crucial Equation.
• SGU 106.
• Codeforces - Ebony and Ivory.
• Codechef - Get AC in one go.
• LightOj - Solutions to an equation.

IV. Tài liệu tham khảo

• https://cp-algorithms.com/algebra/linear-diophantine-equation.html
• Competitive Programming 3 - Steven Halim & Felix Halim.
• https://mathworld.wolfram.com/DiophantineEquation.html