Розширений алгоритм Евкліда

Якщо алгоритм Евкліда обчислює лише найбільший спільний дільник (НСД) двох невід'ємних цілих чисел $a$ і $b$, то розширена версія додатково знаходить спосіб подати НСД через $a$ і $b$, тобто коефіцієнти $x$ і $y$, для яких:

$$a \cdot x + b \cdot y = \gcd(a, b)$$

Важливо зауважити, що за тотожністю Безу таке подання завжди можна знайти. Наприклад, $\gcd(55, 80) = 5$, тому ми можемо подати $5$ як лінійну комбінацію доданків $55$ і $80$: $55 \cdot 3 + 80 \cdot (-2) = 5$

Загальнішу форму цієї задачі розглянуто у статті про лінійні діофантові рівняння. Вона ґрунтуватиметься на цьому алгоритмі.

Коли підходить цей алгоритм?

• Чи окрім $\gcd(a,b)$ вам потрібні самі коефіцієнти $x, y$ такі, що $ax + by = \gcd(a,b)$? (якщо потрібен лише НСД → Алгоритм Евкліда)
• Чи кінцева мета — обернений елемент за модулем $a^{-1} \bmod m$ для необов'язково простого $m$? (інакше для простого $m$ простіше через мала теорема Ферма)
• Чи ви розв'язуєте лінійне діофантове рівняння $ax + by = c$? (→ Лінійні діофантові рівняння)

Алгоритм

У цьому розділі позначатимемо НСД чисел $a$ і $b$ через $g$.

Зміни до вихідного алгоритму дуже прості. Якщо пригадати алгоритм, то можна побачити, що він завершується при $b = 0$ і $a = g$. Для цих параметрів коефіцієнти знаходяться легко, а саме $g \cdot 1 + 0 \cdot 0 = g$.

Починаючи з цих коефіцієнтів $(x, y) = (1, 0)$, ми можемо рухатися назад вгору рекурсивними викликами. Усе, що нам потрібно, — з'ясувати, як змінюються коефіцієнти $x$ і $y$ під час переходу від $(a, b)$ до $(b, a \bmod b)$.

Припустимо, що ми знайшли коефіцієнти $(x_1, y_1)$ для $(b, a \bmod b)$:

$$b \cdot x_1 + (a \bmod b) \cdot y_1 = g$$

і хочемо знайти пару $(x, y)$ для $(a, b)$:

$$a \cdot x + b \cdot y = g$$

Ми можемо подати $a \bmod b$ так:

$$a \bmod b = a - \left\lfloor \frac{a}{b} \right\rfloor \cdot b$$

Підставляючи цей вираз у рівняння коефіцієнтів для $(x_1, y_1)$, отримуємо:

$$g = b \cdot x_1 + (a \bmod b) \cdot y_1 = b \cdot x_1 + \left(a - \left\lfloor \frac{a}{b} \right\rfloor \cdot b \right) \cdot y_1$$

а після перегрупування доданків:

$$g = a \cdot y_1 + b \cdot \left( x_1 - y_1 \cdot \left\lfloor \frac{a}{b} \right\rfloor \right)$$

Ми знайшли значення $x$ і $y$:

$$\begin{cases} x = y_1 \\ y = x_1 - y_1 \cdot \left\lfloor \frac{a}{b} \right\rfloor \end{cases}$$

Реалізація

C++

int gcd(int a, int b, int& x, int& y) {
    if (b == 0) {
        x = 1;
        y = 0;
        return a;
    }
    int x1, y1;
    int d = gcd(b, a % b, x1, y1);
    x = y1;
    y = x1 - y1 * (a / b);
    return d;
}

Python

def gcd(a: int, b: int) -> tuple[int, int, int]:
    # Повертає (g, x, y), де a*x + b*y = g
    if b == 0:
        return a, 1, 0
    d, x1, y1 = gcd(b, a % b)
    # Перераховуємо коефіцієнти під час повернення вгору
    return d, y1, x1 - y1 * (a // b)

TypeScript

// Повертає [g, x, y], де a*x + b*y = g
function gcd(a: number, b: number): [number, number, number] {
  if (b === 0) {
    return [a, 1, 0];
  }
  const [d, x1, y1] = gcd(b, a % b);
  // Перераховуємо коефіцієнти під час повернення вгору
  return [d, y1, x1 - y1 * Math.floor(a / b)];
}

Go

// gcd повертає (g, x, y), де a*x + b*y = g
func gcd(a, b int) (int, int, int) {
    if b == 0 {
        return a, 1, 0
    }
    d, x1, y1 := gcd(b, a%b)
    // Перераховуємо коефіцієнти під час повернення вгору
    return d, y1, x1 - y1*(a/b)
}

Наведена вище рекурсивна функція повертає НСД, а значення коефіцієнтів записує у x і y (які передаються у функцію за посиланням).

Ця реалізація розширеного алгоритму Евкліда дає правильні результати і для від'ємних цілих чисел.

Ітеративна версія

Розширений алгоритм Евкліда також можна записати ітеративно. Оскільки в такому варіанті немає рекурсії, код працюватиме трохи швидше за рекурсивний.

C++

int gcd(int a, int b, int& x, int& y) {
    x = 1, y = 0;
    int x1 = 0, y1 = 1, a1 = a, b1 = b;
    while (b1) {
        int q = a1 / b1;
        tie(x, x1) = make_tuple(x1, x - q * x1);
        tie(y, y1) = make_tuple(y1, y - q * y1);
        tie(a1, b1) = make_tuple(b1, a1 - q * b1);
    }
    return a1;
}

Python

def gcd(a: int, b: int) -> tuple[int, int, int]:
    # Повертає (g, x, y), де a*x + b*y = g
    x, y = 1, 0
    x1, y1, a1, b1 = 0, 1, a, b
    while b1:
        q = a1 // b1
        x, x1 = x1, x - q * x1
        y, y1 = y1, y - q * y1
        a1, b1 = b1, a1 - q * b1
    return a1, x, y

TypeScript

// Повертає [g, x, y], де a*x + b*y = g
function gcd(a: number, b: number): [number, number, number] {
  let x = 1;
  let y = 0;
  let x1 = 0;
  let y1 = 1;
  let a1 = a;
  let b1 = b;
  while (b1 !== 0) {
    const q = Math.floor(a1 / b1);
    [x, x1] = [x1, x - q * x1];
    [y, y1] = [y1, y - q * y1];
    [a1, b1] = [b1, a1 - q * b1];
  }
  return [a1, x, y];
}

Go

// gcd повертає (g, x, y), де a*x + b*y = g
func gcd(a, b int) (int, int, int) {
    x, y := 1, 0
    x1, y1, a1, b1 := 0, 1, a, b
    for b1 != 0 {
        q := a1 / b1
        x, x1 = x1, x-q*x1
        y, y1 = y1, y-q*y1
        a1, b1 = b1, a1-q*b1
    }
    return a1, x, y
}

Якщо уважно придивитися до змінних a1 і b1, то можна помітити, що вони набувають точнісінько таких самих значень, як і в ітеративній версії звичайного алгоритму Евкліда. Отже, алгоритм принаймні обчислить правильний НСД.

Щоб зрозуміти, чому алгоритм обчислює правильні коефіцієнти, зауважимо, що в будь-який момент (перед початком циклу while і наприкінці кожної ітерації) виконуються такі інваріанти:

$$x \cdot a + y \cdot b = a_1$$ $$x_1 \cdot a + y_1 \cdot b = b_1$$

Позначмо значення наприкінці ітерації штрихом ($'$) і припустимо, що $q = \frac{a_1}{b_1}$. З алгоритму Евкліда маємо:

$$a_1' = b_1$$ $$b_1' = a_1 - q \cdot b_1$$

Щоб виконувався перший інваріант, має бути правильно, що:

$$x' \cdot a + y' \cdot b = a_1' = b_1$$ $$x' \cdot a + y' \cdot b = x_1 \cdot a + y_1 \cdot b$$

Аналогічно для другого інваріанта має виконуватися:

$$x_1' \cdot a + y_1' \cdot b = a_1 - q \cdot b_1$$ $$x_1' \cdot a + y_1' \cdot b = (x - q \cdot x_1) \cdot a + (y - q \cdot y_1) \cdot b$$

Порівнюючи коефіцієнти при $a$ і $b$, можна вивести рівняння оновлення для кожної змінної, які гарантують, що інваріанти зберігаються протягом усього алгоритму.

Наприкінці ми знаємо, що $a_1$ містить НСД, тож $x \cdot a + y \cdot b = g$. Це означає, що ми знайшли потрібні коефіцієнти.

Код можна оптимізувати ще більше: прибрати змінні $a_1$ і $b_1$ та просто повторно використати $a$ і $b$. Однак якщо так зробити, то ми втратимо можливість міркувати про інваріанти.

Задачі для практики

• UVA - 10104 - Euclid Problem
• GYM - (J) Once Upon A Time
• UVA - 12775 - Gift Dilemma

Відеоматеріали