밀러-라빈 소수 판별법

소수 판별법

정보 문제에서 심심치 않게 등장하는 것이 바로 소수 판별법이다. 또, 소수 판벌법은 그 자체로도 중요하지만, 다른 정수론 문제에서 기본이 되는 만큼 그 효율이 중요하다.

흔히 생각할 수 있는 방법은 2에서 n-1까지의 수로 n을 나누어 보는 것이다. 이 경우 n이 커지면 시간이 너무 많이 걸린다. 조금 더 효율적으로는 2에서 \(\sqrt{n}\)까지의 수로 나누는 것이다. \(\sqrt{n}\)보다 큰 n의 약수 m이 있다면 이미 \(n/m\)에서 확인했기 때문이다.

또 다른 방법으로는 에라토스테네스의 체이다. 이 방법은 n까지의 소수를 전부 찾을 수 있다는 장점이 있다. 그러나 아직 느리다. 더 빠른 알고리즘이 필요하다.

밀러-라빈 소수 판별법

더 빠른 알고리즘이 여기 있다. 바로 밀러-라빈 소수 판별법이다. 문제는 이 알고리즘이 정확하지 않다는 것이다. 무려 오류가 있다! 그러나 그 오류는 수가 엄청 클 때 발생한다. 정확히는 우리가 흔히 사용하게 된 '작은 수', 특히 C기준 long long에서는 그 예외가 없다는 것이다. 그래서 웬만한 문제에서는 이를 사용해도 된다.

왜 오류가 있을까? 바로 밀러-라빈은 소수이기 위한 충분 조건을 찾아주는 알고리즘이 아니다. 필요 조건을 찾아준다. 그러나 그 예외가 아주 큰 수에서 처음 등장하는 것뿐이다.

이번 글에서는 밀러-라빈의 이론적 배경과 알고리즘의 구현에 대해 알아보자.

밀러-라빈 소수 판별법의 이론적 배경

우선 밀러-라빈을 이해하기 위해서는 정수론에 대한 기본적인 지식이 필요하다. 어려울 것 같다고 겁먹지 말자. 우리는 수학을 하는 것이 아니므로 수학적으로 어려운 내용을 다루지는 않을 것이다. 딱 밀러-라빈에 필요한 내용만 알짜배기로 배우자.(사실 필자도 잘 모른다. 수전 해줘.)

Lemma 1 - 페르마의 소정리

그 유명한 페르마의 마지막 정리...가 아닌 페르마의 소정리다.

\(a^{p}\equiv a\pmod p\)

\(a^{p-1}\equiv 1\pmod p\)

이때 a는 정수, p는 소수다. 복잡하게 생각하지 말자. 유명한 수학자들이 증명해주었으니 우리는 이용만 하자.

Lemma 2

\(x^2 \equiv 1\pmod p \Rightarrow x \equiv 1 \pmod p \ or \ x \equiv -1 \pmod p\)

이건 어렵지 않게 증명할 수 있으나 그냥 넘어가도록 하자. (궁금하면 댓글로, 아님 직접 증명...)

밀러-라빈 소수 판별법

위의 두 Lemma를 이용해서 밀러-라빈의 이론적 배경을 알아보자. n이 소수인지 판별해보자.

For all prime number n, if \(n>2\), then

\(n-1=2^{t}m, \ ( \text{t is an nonpositive integer , m is an odd integer}) \)

From Lemma 1, for positive integer a,

\(a^{2^t m} \equiv 1\pmod n \ (\because a^{n-1} \equiv 1 \pmod n)\)

Thus,

\(a^{(2^{t-1} m)\times 2} \equiv 1\pmod n\)

From Lemma 2,

\(a^{2^{t-1} m} \equiv 1 \pmod n \text{ or }a^{2^{t-1} m} \equiv -1 \pmod n\)

In the first case, we can continually use Lemma 2 and get the following result.

\(\exist k<t, \ a^{2^{k} m} \equiv -1 \pmod n, \text{ k is an positive integer}\)

or

\(a^{m} \equiv 1 \pmod n \text{ or } a^{m} \equiv -1 \pmod n \)

정리

따라서 n이 소수라면 다음이 성립한다.

또, long long 범위에서는 다음 수에 대해서만 성립하면 소수라는 것이 증명되었다.

밀러-라빈 소수 판별법의 구현

def power(a,k,n):
    if k==1:
        return a%n
    else:
        tmp=power(a,k//2,n)
        if k%2==0:
            return (tmp*tmp)%n
        else:
            return (tmp*tmp*a)%n
checklist=[2,3,5,7,11,13,17,19,23,29,31,37]
def isprime(n):
    if n==2 or n==3:
        return 1
    if n%2==0:
        return 0
    c=1
    for a in checklist:
        if a%n==0:
            return 1
        d=n-1
        while True:
            t=power(a,d,n)%n
            if t==n-1:
                break
            if d%2==1:
                if t==1 or t==n-1:
                    break
                else:
                    c=0
                    break
            d=d//2
    return c