GCD & Extended Euclidean Algorithm

Greatest Common Divisor(GCD)

두 수를 모두 나누는 가장 큰 수

최대공약수가 무엇인지는 모두 알 것이다.
PS에서도 생각보다 자주 쓰이는 개념이다.

• • • Greatest Common Divsior
    • Euclidean Algorithm
    • Extended Euclidean Algorithm

Greatest Common Divsior

먼저, 정의를 조금 다르게 해보자.
$$gcd(a,b) = \begin{cases} 0 & (a=0, b=0) \\ b & (a=0, b > 0) \\ a & (a > 0, b = 0) \\ \max \limits_{d \vert a, d \vert b} d & (a > 0, b > 0) \end{cases}$$

둘 다 양수일 때는 원래의 정의와 동일하다.
둘 중 하나가 0이고, 하나가 0이 아닐 때도 잘 생각해보면 원래의 정의에 부합한다는 것을 알 수 있다.
하지만 $gcd(0,0) = 0$이라는 것은 원래의 정의에 어긋난다. 0이 0을 나누고, 0보다 큰 수는 0을 나누지 않는다는 말이 되기 때문이다.
하지만 이렇게 정의함으로서 얻는 장점들이 있다.

예를 들어, 어떤 배열의 모든 원소의 gcd를 구하는 코드를 이렇게 작성할 수 있다.

int g = 0;
for(int i = 0; i < n; ++i) g = gcd(g, A[i]);
cout << g << endl;

만약 A 배열이 {6, 3, 15}라면 g=3이 출력될 것이다.
그런데 A 배열이 {0, 0, 6, 0, 3, 15, 0}이라면?
gcd(0, 0) = 0이고, gcd(0, a) = a이기 때문에, 0들은 무시되고 {6, 3, 15}의 gcd를 구한 것과 같은 결과가 나온다.
마치 덧셈에서의 0 혹은 곱셈에서의 1과 같은, 항등원으로서의 역할을 할 수 있는 것이다.
이는 각종 자료구조를 초기화하거나 임시 변수를 초기화할 때의 문제를 해결해준다.

Euclidean Algorithm

유클리드 알고리즘은 두 수의 최대공약수를 구하는 알고리즘으로, 가장 오래된 알고리즘 중 하나라고도 알려져있다.
이 알고리즘이 사용하는 핵심 원리는 다음과 같다.

$b > 0$이라고 하자.
산술의 기본 정리에 의해 $a=bq+r \quad (0 \leq r < b)$인 $r$이 존재한다.
이 때, $gcd(a, b) = gcd(b, r)$.

이에 대한 증명은 어렵지 않으니 생략하도록 한다.
또한 위의 정의에 따르면 $b=0$일 때 $gcd(a,b)=a$이다.
따라서…

$$gcd(a,b) = \begin{cases} gcd(b, r) \quad & (0 \leq r < b, a=bq+r) \\ a \quad & (b=0) \end{cases}$$

단 한 줄로 코딩할 수 있다.

int gcd(int a, int b){ return b?gcd(b,a%b):a; }

재귀 함수라서 시간이 걱정된다면 그럴 필요 없다. tail recursion이 일어나는 대표적인 경우이기 때문에, 최적화 옵션이 켜져있다면 반복문으로 처리된다! (재귀함수 코드와, -O2 옵션부터 tail recursion에 의해 반복문으로 처리된 코드의 컴파일 결과를 보자. 반복문을 직접 짠 코드와 어셈블리가 동일하다!)

Extended Euclidean Algorithm

먼저 둘 중 하나는 0이 아닌, 음이 아닌 정수 $a$와 $b$를 생각하자.
$ax+by=gcd(a,b)$가 되는 정수 $(x, y)$를 찾는 문제가 있다. (이렇게 정수 해만을 찾아야 한다는 조건이 붙은 방정식을 Diophantine equation이라고 부른다.)
한 해(특수해 라고 부르기도 한다)를 알고 있으면 이 방정식의 모든 해를 표현하는 식을 만들 수 있다. (어떤 해 $(x,y)$를 찾고 나면, 당장 $(x+b, y-a)$ 또한 해가 될 수 있다! 일반해는 이보다 약간 더 복잡하다.)
하지만 한 해를 찾는 것 자체가 쉽지 않다. $0 \leq i < b$ 사이의 모든 $i$만 확인해도 되긴 하지만, $b$가 너무 커서 확인할 수 없을 수도 있다.
그럴 때 유클리드 알고리즘과 같은 원리로 특수해를 하나 찾을 수 있다. 생각해보자.
위의 코드처럼 재귀적으로 뭔가를 구하고, 그 결과를 이용해서 자신의 답을 구한다고 하자.
우리는 $ax+by = gcd(a,b)$가 되는 $(x,y)$를 찾고 싶다.
그런데 우리가 알고 있는 것은 $a=bq+r$일 때, $bx' + ry' = gcd(b,r) = gcd(a,b)$가 되는 $(x', y')$이다.
$bx'+ry'$의 식에서 $r$을 없애면 $(x,y)$를 찾을 수 있지 않을까? 그를 위해 식을 변형해보자.
$r=a-bq$ 이므로 $bx'+(a-bq)y'=gcd(b,r)=gcd(a,b)$ 이다.
즉 $a \cdot y' + b \cdot (x' - qy') = gcd(a,b)$ 이다.
또한, $b=0$이면 $a=gcd(a,b)$이므로 $a \cdot 1 + b \cdot 0 = gcd(a,b)$이다. 즉 이 때의 해는 $(1,0)$이다.

std::pair를 이용하여 다음과 같은 코드를 작성할 수 있다.

using namespace std;
pair<int,int> ext_gcd(int a,int b){
    if(b){
        auto tmp = ext_gcd(b, a%b);
        return {tmp.second, tmp.first - (a/b) * tmp.second};
    } else return {1, 0};
}

확장 유클리드 알고리즘이 필요한 대표적인 경우로 modular inverse가 있다.
어떤 $a$와 $M$이 주어질 때, $ax \equiv 1 \mod M$인 $x$를 찾는 문제이다.
이러한 $x$는 $gcd(a, M)=1$일 때만 존재하기 때문에, 대부분의 문제에서는 $M$이 소수인 경우에 사용하게 된다.
저런 $x$를 위의 Diophantine equation 형태로 바꾸려면 어떻게 해야 할까?
$ax \equiv 1 \mod M$ 이라는 식은 $ax + My = 1$인 정수 $y$가 존재한다는 뜻이다. 즉 저 식의 어떤 해 $x$를 찾으면 된다.
음수가 나오는 것을 막기 위해 이렇게 코딩할 수 있다.

int mod_inv(int a, int M){
    return (ext_gcd(a, M).first + M) % M;
}

modular inverse를 제대로 짰는지 확인하는 것은 매우 쉽다.
랜덤한 숫자 x를 잡아보며, x * mod_inv(x, M) % M == 1 인지 확인하면 된다.

modular inverse의 다른 방법으로는 페르마 소정리(Fermat’s little theorem)를 이용한 것이 있다.
$a^{ \phi (M) } \equiv 1 \mod M$ 이므로, $M$이 소수이면 $a \cdot a^{M-2} \equiv 1$이다. 즉 $a^{M-2}$를 modular inverse로 사용할 수 있다.

Written with StackEdit.