본문 바로가기

알고리즘/수학

소수 판별법 - 에라토스테네스의 체, 밀러-라빈(Miller-Rabin) 소수판별법

728x90
반응형

어떠한 자연수 N이 소수인지를 판별하는 방법은 여러 가지 방법이 있다. 

그중에서 너무 난도 높은 것은 제외하고 충분히 PS에서 쓸만한 방법을 알아보자. 

우선, N이 소수인지를 판별하는 경우N이하의 소수가 몇개있는지, N이하의 소수를 모두 구하는 경우 두가지로 보통 나뉜다. 

 

  1. N을 2부터 N-1까지 모두 나누기 

소수는 1또는 자기 자신으로만 나누어지기 때문에 위의 방법으로 N이 나누어 떨어지지 않는다면 N은 소수라고 할 수 있다. 

만약, N이 어떤 두 수로 나누어진다고 한다면, N = p*q의 형태가 되고, p와 q 둘 중 하나는 반드시 $ \sqrt{N} $의 값 이하일 것이다. 

따라서 N-1까지 나눌 필요는 없고, $ \sqrt{N} $까지만 나눠도 소수인지 판별이 가능하다. 

물론, 이 방법을 이용하도록 출제되는 문제는 거의 없을 것이다.... 시간 복잡도가 $ O(\sqrt{N}) $이기 때문에 N이 커지면 소수를 판별하는데 시간이 오래 걸린다.

또, N이하의 소수를 모두 구하라고 한다면 $ O(N\sqrt{N}) $의 시간복잡도를 갖기 때문에 N의 범위가 매우 제한적이다. 

1
2
3
4
5
6
7
8
9
10
int main(void){
    int N;
    for(int i = 2; i*i<N; i++){
        if(N%i == 0) {
            cout<<"N is not prime";
            return 0;
        }
    }
    cout<<"N is prime";
}

 

  2. 에라토스테네스의 체

소수의 개수 또는 소수들을 구하는 방법 중에서 가장 보편적으로 이용되는 방법이다. 

초기값이 0인 N개의 원소를 갖는 배열(A)을 선언하여, 소수가 아닌 수는 체크해주는 방식이다.  

방식은, i = 2부터 i = N까지 순서대로 돌면서, 만약 A[i]가 0이라면 A[k*i] (k = 2, 3, 4, ...) 들은 모두 1을 할당해준다. 

그렇게 되면 결국 i가 소수인 경우에만 A[i] = 0이 되고, 나머지는 모두 1이 되어버린다. 

소수의 배수들은 모두 합성수라는 성질을 이용하여 "체로 걸러내듯이" 합성수를 모두 걸러내는 방식이다. 

시간 복잡도는 $ O(nlog(logn)) $ 으로, 보통 문제를 풀 땐 100만까지 정도 소수를 구해낸다. 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
#include<iostream>
using namespace std;
int main(void) {
    int arr[1000= { 0, };
    for (int i = 2; i < 1000; i++) {
        if (arr[i] == 0) {
            for (int j = 2; j*< 1000; j++) {
                arr[i*j] = 1;
            }
        }
    }
    for (int i = 2; i < 1000; i++) {
        if (arr[i] == 0cout << i << ' ';
    }
}
 

 

코드를 통해서도 알 수 있듯이, 이 방식은 어떤 특정한 수가 소수인지를 판단하기에는 좋은 방법이 아니다.

N이하의 소수가 몇개 있는지, 소수들을 직접 구하는 경우에 유용하게 쓰일 수 있다. 

 

  3. 밀러-라빈 소수 판별법 (Miller-Rabin Primality Test)

밀러-라빈 소수 판별법은 특정한 수 N이 소수인지 '확률적으로' 판단할 수 있는 판별법이다.

'확률적으로' 라는 의미는, 밀러-라빈 소수 판별법을 이용했을 때, N이 '합성수'이거나, '아마도 소수일 것이다'라고 판단을 해준다는 것이다.  

대신 우리가 생각하는 큰 수 (long long) 범위 내에서는 소수인지 정확히 판단이 가능하므로 너무 걱정할 필요는 없을 것 같다.

 

밀러 - 라빈 소수판별법을 이해하기 위해서는 정수론 개념이 필요하다. 

우선 이 판별법을 위한 보조 정리가 이용된다. 

 

[보조 정리 1] 페르마의 소정리 (Fermat's Little Theorem)

소수 $p$와 정수 $a$에 대해서 다음을 만족한다. 

$a^p \equiv a\:(mod\: p) $

$a^{p-1} \equiv 1\:(mod\: p)\quad if \quad a\%p \neq 0 $

 

페르마의 소정리에 대한 증명은 이 사이트를 참고하면 도움이 될 것 같다. https://casterian.net/archives/382

 

페르마의 소정리 (Fermat’s Little Theorem) – The Casterian

페르마의 소정리) 소수 $p$와 정수 $a$에 대해 \[ a^p \equiv a \pmod p \] 아주 간단하지만 강력한 정리입니다. 보통 문제에서 어떤 수의 거듭제곱을 다른 어떤 수로 나눈 나머지를 구하는 과정이 있는데

casterian.net

 

[보조 정리 2]

소수 $p$에 대해서 $ x^2 \equiv 1 \: (mod \: p) $ 이면 , $ x \equiv -1 \: (mod \: p) \quad or \quad x \equiv 1 \: (mod \: p) $ 이다.   

 

( $x^2 - 1 $이 $p$의 배수이므로, $x-1$ 또는 $x+1$ 이 $p$의 배수가 된다. 따라서 위의 보조정리가 성립한다. )

이 두 보조정리를 이용하여 소수를 판별할 수 있다. 

$N$이 2가 아닌 소수라고 가정하자. 그러면 $N-1$은 짝수가 되고, $N-1$을 홀수가 될 때까지 2로 나눠주면 다음과 같이 표현이 가능하다. 

$N-1 = d\times 2^r $

$N$이 소수이기 때문에 양의 정수 $a$에 대해서 페르마의 소정리를 만족하게 된다. 따라서 다음 식을 만족한다.

$a^{N-1} \equiv 1 \: (mod \:N) $ , $a^{d\times 2^r} \equiv 1 \: (mod \:N) $ 

즉, $(a^{d\times 2^{r-1}})^2 \equiv 1 \: (mod \:N) $

따라서 [보조정리 2]에 의해서 $a^{d\times 2^{r-1}} \equiv -1 \: (mod \:N) \quad or \quad a^{d\times 2^{r-1}} \equiv 1 \: (mod \:N) $ 를 만족한다. 

만약 후자인 경우에는 또다시 [보조 정리 2]를 이용할 수 있다. 

계속해서 1이 나오는 경우에는 결국 $ a^{d\times 2} \equiv 1 \: (mod \:N) $ 까지 나와서 마지막으로

$ a^d \equiv 1 \: (mod \: N) \quad or \quad a^d \equiv -1 \: (mod \: N) $ 의 결과가 나오게 된다.

여기서는 d가 홀수이기 때문에 더이상 [보조 정리 2]를 사용하지 못한다. 

위 과정을 통해서 $N$이 소수라면 $N$보다 작은 양의 정수 $a$에 대해서 다음 둘 중 하나를 만족한다고 할 수 있다. 

1. $a^d \equiv 1 \: (mod \: N) $

2. $a^{d \times 2^r} \equiv -1 \: (mod \: N) \quad for\; some \;r \quad (0 \leq r < s) $

 

여기서 $N$이 '아마도 소수일 것이다'라고 확률적으로 말하는 이유는 $N$이 소수가 아니더라도 특정한 a에 대해서 이를 만족할 수 있기 때문이다.

따라서 $N$이 '소수이다' 라고 단정 짓기 위해서는 최대한 많은 a를 적용시켜보아야 한다. 

int 범위의 $N$을 판별하기 위해서는 a = 2, 7, 61 세 수만 통과하면 $N$을 소수라고 할 수 있고, long long 범위의 $N$을 판별하기 위해서는 37까지의 소수들에 대해서만 통과하면 $N$을 소수라고 할 수 있다고 알려져 있다. 

반대로 $N$이 합성수라면, [보조 정리 1, 2]를 만족하지 못하므로, 중간에 1 또는 -1이 아닌 다른 값이 나머지로 나올 것이다. 

소스 코드)

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
long long power(long long x, long long y, long long mod){
    long long res = 1;
    x %= mod;
    while(y){
        if(y % 2) res = (res * x) % mod;
        y /= 2;
        x = (x*x) % mod;
    }
    return res;    
}
 
//if n is prime, return true
bool miller(long long n, long long a){
    if(a % n == 0return false;
    long long k = n-1;
    while(1){
        long long temp = power(a, k, n);
        if(temp == n-1return true//a^k = -1 (mod n)
        if(k%2return (temp == 1 || temp == n-1);
        k /= 2;
    }
}
 

 

$a^{n-1} $을 계산해주어야 하므로 빠른 거듭제곱 기법(power 함수)을 이용한다. 

또한, overflow를 방지해주기 위해서 계속해서 mod로 나눈 나머지를 이용한다. 

 

예제) BOJ 5615 아파트 임대

문제 : https://www.acmicpc.net/problem/5615

 

5615번: 아파트 임대

문제 동규부동산에서 아파트를 임대하고 있다. 아파트의 방은 아래 그림과 같이 면적이 2xy + x + y이다. (x와 y는 양의 정수) 동규부동산의 카탈로그에는 아파트의 면적이 오름차순으로 적혀져 있�

www.acmicpc.net

이 문제는 결과적으로 어떠한 수가 소수인지를 판별하는 문제이다. 

수의 범위가 int를 넘어가기 때문에 밀러-라빈 소수 판별법을 이용해야 TLE가 뜨지 않는다. 

 

소스코드)

더보기
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
#include<iostream>
 
#define ull unsigned long long int
#define FOR(i, N) for(int i = 1; i<=N; i++)
 
using namespace std;
ull miller_base[] = { 2,761 };
ull power(ull x, ull y, ull mod);
bool miller(ull n, ull a);
 
int main(void) {
    ios::sync_with_stdio(false); cin.tie(NULL); cout.tie(NULL);
 
    int N; cin >> N;
    ull a;
    int ans = 0;
    FOR(j, N){
        cin >> a;
        a = a * 2 + 1;
        FOR(i, 3) {
            if (miller(a, miller_base[i - 1]) == false) {
                ans++;
                break;
            }
        }
    }
    cout << N- ans;
}
ull power(ull x, ull y, ull mod) {
    ull res = 1;
    x %= mod;
    while (y) {
        if (y % 2) res = (res*x) % mod;
        y /= 2;
        x = (x*x) % mod;
    }
    return res;
}
 
bool miller(ull n, ull a) {
    if (a % n == 0return true;
    ull k = n - 1;
    while (1) {
        ull temp = power(a, k, n);
        if (temp == n - 1return true
        if (k % 2return (temp == 1 || temp == n-1);
        k /= 2;
    }
}
 

 

 

PC로 보시는 것을 권장합니다. 

피드백은 언제나 환영입니다. 댓글로 달아주세요 ^-^

728x90
반응형
  • 이재현 2021.05.16 22:28 댓글주소 수정/삭제 댓글쓰기

    수의 범위가 int를 넘어가는데 2, 7, 61만 봐도 되나요?

    • 아파트 임대 문제를 말씀하시는거라면, 정확히는 int보다 조금 더 큰 unsigned int의 범위, 즉 약 2^32 정도까지 2, 7, 61로 확인가능합니다.
      자세한건 https://ko.m.wikipedia.org/wiki/%EB%B0%80%EB%9F%AC-%EB%9D%BC%EB%B9%88_%EC%86%8C%EC%88%98%ED%8C%90%EB%B3%84%EB%B2%95 를 참고해주세요 :)