소인수분해/알고리즘

  상위 문서: 소인수분해

1. 개요[편집]

---

암호학이 대두되면서 큰 소수를 구할 필요가 생겼고, 그와 동시에 큰 합성수의 소인수분해에 대한 연구도 진행되었다.

2. 알고리즘[편집]

---

2.1. 기초[편집]

---

주어진 숫자 n의 소인수를 구한다고 할 경우 아래와 같은 순서로 진행하면 된다.
어차피 작은 소수에서 n%i == 0이 성립하지 못한다면 그보다 큰 합성수는 n%i == 0가 성립하지 못하므로 n%i == 0이 성립하는 첫번째 i는 소수라는 점을 이용한 알고리즘이다.
이 알고리즘으로 소인수를 구하면 천억이 넘는 숫자도 소인수가 순식간에 구해진다.

간단하게 파이썬으로 코딩 해보면 다음과 같다.

# -*- coding: utf-8 -*-
import os
import sys

def find_prime(input_num):
    if input_num <= 2:
        return [input_num,]
    for idx in range(2,input_num):
        if input_num % idx == 0:
            ret_list = []
            val_a    = find_prime(idx)
            val_b    = find_prime(int(input_num/idx))
            ret_list = val_a + val_b
            return ret_list
    return [input_num,]

def main():
    try:
        input_num  = int(sys.argv[1])
    except:
        print("usage: python main.py <number>")
        print(" ex>   python main.py 12345")
        quit()
    prime_list = find_prime(input_num)
    check_num  = 1
    for prime_at in prime_list:
        check_num *= prime_at
    print("[%s] input_num=%d, check_num=%d" % (input_num == check_num, input_num, check_num))
    print("%s" % prime_list)

main()

다만 위의 프로그램은 프로그래밍을 익히는 용도로는 유용할수 있으나, 실제 소인수분해 용으로 쓰기에는 적합하지 않다. 실제로 컴퓨터로 다루는 수는 1000억은 '겨우'라는 소리가 나올만큼 큰 수를 다룰 필요가 있기 때문이다.

가장 쉬운 방법은 다른 프로그램을 이용해서 미리 소수 테이블을 작성해 두고, 이를 활용하는 것이다. 예를 들어 2¹⁶ 보다 작은 소수는 6542개인데, 이를 미리 배열에 저장해 두면, 42억 (=2³²) 보다 작은 수는 겨우 6542번 나누어 보기만 하면 소수인지 판정하거나, 1개 이상의 소인수를 구할 수 있다. 2³²보다 작은 소수는 모두 약 2억개(203,280,221개)인데, 이를 미리 구해서 적절한 DB에 저장해 두면 약 1844경 (= 2⁶⁴ = 18,446,744,073,709,551,616)보다 작은 수의 소인수 분해를 쉽게 할 수 있다.

20자리 정도 되는 이정도 수면 큰 수라고 생각할 수 있지만, RSA 암호화같은 암호학에서 기본 수백자리 수를 다뤄야 한다. RSA 넘버를 보면 가장 작은 것이 100자리 부터 시작하는데, 그마저 250자리 수 까지는 모두 소인수분해되었다. 가장 큰 RSA-2048은 무려 617자리 수이다. 수의 단위 중 이름이 붙은 최고 단위인 구골이 101자리 수라는 것을 생각해보면, RSA 암호화에서 다루는 수는 사람에게는 어마어마하게 큰 수이다.

밑과 지수를 구하는 기본적인 파이썬 함수는 여기서 설명되어 있다.

2.2. 페르마 소인수분해법[편집]

---

곱셈 공식을 응용한 소인수분해법.
[math(N=x^2-y^2)]을 만족하는 두 자연수 [math(x)], [math(y)]를 찾을 수 있다면, [math(x+y)]와 [math(x-y)] 모두 [math(N)]의 약수가 된다. 식을 조금 변형하여, [math(N+y^2=x^2)]이 제곱수가 되도록 하는 자연수 [math(y)]를 찾는 문제가 된다.
특성상 [math(\sqrt{N})]에 가까운 약수를 찾을 때 유용하며, [math(N)]이 5 이상의 홀수일 때만 적용 가능하다. 또한, [math(N)]이 소수일 경우 성능이 매우 떨어지기 때문에, [math(y)]의 값이 어느 정도 커지면 위의 기초 알고리즘을 대신 사용해야 한다.
한두 개의 소인수만 발견할 수 있기 때문에, 소인수분해를 마치려면 아래 알고리즘을 재귀적으로 수행해야 한다.
유사 코드로 표현하면 다음과 같다. 변수 [math(x_2)], [math(y_2)]는 곱셈 연산을 한 번으로 줄이기 위해 도입한 변수이다.

1. 홀수 [math(N)]을 입력받는다.
2. [math(x=\left\lceil\sqrt{n}\right\rceil)], [math(y=0)]
3. [math(x_2=x^2)], [math(y_2=0)]
4. [math(N+y_2)]의 값을 [math(x_2)]와 비교하여, 같으면 5로, 좌변이 크면 6로, 우변이 크면 8로.
5. [math(x-y)]를 출력하고 종료한다. 이는 [math(N)]의 비자명한 약수이다.
6. [math(x_2\leftarrow x_2+2x+1)], [math(x\leftarrow x+1)]
7. 5로
8. [math(y_2\leftarrow y_2+2y+1)], [math(y\leftarrow y+1)]
9. [math(y\geq N/6)]일 때, [math(N)]은 소수이다. 알고리즘을 종료한다.
10. 5로

[math(N)]을 4로 나눈 나머지가 1일 때, [math(x)]는 홀수, [math(y)]는 짝수이고, [math(N)]을 4로 나눈 나머지가 3일 때, [math(x)]는 짝수, [math(y)]는 홀수이다. 이를 이용하여 실행 시간을 절반으로 줄일 수 있다. 기초 알고리즘과 결합하면, 다음과 같다.

1. 홀수 [math(N)]을 입력받는다.
2. [math(x=\left\lceil\sqrt{n}\right\rceil)], [math(y=y_2=0)]
3. [math(N\equiv3\ (\mathrm{mod}\ 4))]일 때, [math(y=y_2=1)]
4. [math(x+y)]가 짝수일 때, [math(x\leftarrow x+1)]
5. [math(x_2=x^2)]
6. [math(N+y_2)]의 값을 [math(x_2)]와 비교하여, 같으면 7로, 좌변이 크면 8로, 우변이 크면 10으로.
7. [math(x-y)]를 출력하고 종료한다. 이는 [math(N)]의 비자명한 약수이다.
8. [math(x_2\leftarrow x_2+4x+4)], [math(x\leftarrow x+2)]
9. 7로
10. [math(y_2\leftarrow y_2+4y+4)], [math(y\leftarrow y+2)]
11. [math(2y<x)]일 때, 8로. 아래 코드는 [math(\sqrt{N}/3)] 이하의 소인수를 찾기 위한 기초 소인수분해 알고리즘이다.
12. [math(i=3)], [math(i_2=27)]
13. [math(N)]이 [math(i)]의 배수일 때, [math(i)]를 출력하고 종료한다.
14. [math(i_2\leftarrow i_2+12i+12)], [math(i\leftarrow i+2)]
15. [math(i_2>N)]일 때, [math(N)]은 소수이다. 알고리즘을 종료한다.
16. 13으로

2.3. 오일러 소인수분해법[편집]

---

2.4. 폴라드 로 알고리즘[편집]

---

1975년 존 폴라드가 발표한 소인수분해 알고리즘. 확률적 알고리즘으로, [math(N)]이 합성수더라도 소인수분해에 가끔 실패할 수 있다. 다만 시간복잡도는 확실히 줄어든다.

다음과 같은 함수를 생각하자.

• [math(g(x)=x^2+1\mod\ N)]

[math(g)]는 유한집합 [math(\mathbb{N}_N)]에서 [math(\mathbb{N}_N)]으로의 함수이므로, 다음과 같이 정의된 수열에서 반드시 충돌이 발생한다.

• [math(x_{i+1}=g(x_i))]

[math(N)]의 소인수 [math(p)]를 하나 가정하여, 이번에는 다음과 같은 함수와 수열을 생각하자.

• [math(h(y)=y^2+1\mod\ p)]
• [math(y_0=x_0\mod\ p)], [math(y_{i+1}=h(y_i))]

그러면 [math(y_i=x_i\mod\ p)]가 성립함을 알 수 있으며, 마찬가지로 수열 [math(\left\{y_i\right\})]에서도 반드시 충돌이 발생한다.

[math(x_i=x_j)]이면 [math(y_i=y_j)]이기 때문에, [math(\left\{x_i\right\})]에서 발생한 충돌은 [math(\left\{y_i\right\})]에서도 충돌을 일으킨다. 만약 [math(y_i=y_j)]인데 [math(x_i\neq x_j)]이라면, [math(\gcd(x_i-x_j,N))]은 [math(p)]의 배수이지만 [math(N)]이 아니다. 이러한 [math(i)]와 [math(j)]를 찾는 게 폴라드 로 알고리즘의 핵심이다.

[math(\mathbb{N}_N)]을 정점으로 하고 [math((x,g(x)))]를 간선으로 하는 유향 그래프를 생각하자. 위에서 말한 성질로 인해 이 그래프에는 반드시 사이클이 존재하며, 이는 [math(\left\{x_i\right\})]가 주기성을 띤다는 것을 의미한다. Floyd's cycle-finding algorithm을 이용하여 그 주기를 찾을 수 있다.

[math(x_i=x_{2i})]가 성립하는 가장 작은 자연수 [math(i)]에 대해, 사이클의 주기 [math(d)]는 [math(i)]의 약수이며, 임의의 음이 아닌 정수 [math(j)]에 대해 [math(x_{i+dj}=x_{2(i+dj)})]이다.

[math(\left\{y_i\right\})]의 주기가 [math(\left\{x_i\right\})]의 주기보다 작은지는 수열을 직접 계산하기 전까지는 알 수 없으며, 혹시나 주기가 같다면 [math(x_0)]의 값을 달리하여 알고리즘을 다시 돌리면 된다.

폴라드 로 알고리즘을 유사 코드로 표현하면 다음과 같다.

1. 자연수 [math(N)]을 입력받는다.
2. [math(N=25)]일 때, [math(5)]를 출력한다.[1]
3. [math(x)]를 임의로 정한다.
4. [math(y=x)]
5. [math(x\leftarrow g(x))]
6. [math(y\leftarrow g(g(y)))]
7. [math(x=y)]이면 3로
8. [math(d=\gcd(x-y,N))]
9. [math(d=1)]이면 5로
10. [math(d)]는 [math(N)]의 자명하지 않은 약수이다. [math(d)]를 출력하고 종료한다.

2.5. 이차 체[편집]

---

Quadratic Sieve

페르마 소인수분해법의 단점은 조건을 만족하는 [math(y)]를 찾기 매우 어렵다는 데 있다. 홀수 [math(N)]에 대해 [math(N+y^2=x^2)]을 만족하는 자연수 [math(y)]의 개수는 [math(\sqrt{N})]보다 작은 [math(N)]의 약수의 개수와 같다. 페르마 소인수분해법은 기초 알고리즘에서 나눗셈을 덧셈으로 바꾼 정도밖에 시간복잡도를 개선하지 못했다.
그런데 [math(y)]를 일일이 탐색하는 것이 아니라 만들어낼 수 있다면 어떨까? 이차 체는 이를 실제로 구현한 알고리즘이다.

다음과 같은 함수 [math(Q)]를 생각하자.

• [math(Q(x)=x^2-N)]

[math(Q(x))]가 완전제곱수가 되면 그걸로 끝이지만, 그럴 확률은 너무나도 낮다. 대신 자연수 [math(x_1)], [math(x_2)], ⋯, [math(x_k)]를 몇 개 골라 [math(Q(x_1))], [math(Q(x_2))], ⋯, [math(Q(x_k))]를 계산한 다음, 이 중에서 몇 개를 다시 선별하여 [math(Q(x_{i_1})Q(x_{i_2})\cdots Q(x_{i_r}))]가 완전제곱수가 되도록 할 것이다. 그러면

• [math(Q(x_{i_1})Q(x_{i_2})\cdots Q(x_{i_r})\equiv(x_{i_1}x_{i_2}\cdots x_{i_r})^2\ (\mathrm{mod}\ N))]

이 성립할 것이고, 페르마 소인수분해법에 사용할 [math(y)]를 50% 이상의 확률로 찾을 수 있다.

[math(N)]의 크기에 맞추어 적당한 [math(B)]를 정한다. 모든 소인수가 [math(B)] 이하인 정수를 [math(B)]-smooth number라 한다. 다음과 같은 과정을 거쳐서 수열 [math(\left\{x_i\right\})]을 생성한다.

• [math(n=0,\ \pm1,\ \pm2,\cdots)]에 대해, [math(Q\left(\left\lfloor\sqrt{N}\right\rfloor+n\right))]이 [math(B)]-smooth일 때, [math(\left\lfloor\sqrt{N}\right\rfloor+n)]을 수열 [math(\left\{x_i\right\})]에 추가한다.

수열에 원소를 추가할 때마다, 마지막 원소 [math(x_k)]에 대해 [math(Q(x_k))]를 소인수분해한 결과를 행렬 [math(E)]에 추가한다. 즉, [math(p_0=-1)]이라 하고, [math(B)] 이하의 소수의 집합을 [math(\left\{p_1,\ p_2,\ \cdots,\ p_{\pi(B)}\right\})]라 할 때, [math(\displaystyle Q(x_i)=\prod_{j=0}^{\pi(B)}p_j^{e_{i,\ j}})]가 된다.
[math(E)]에 열을 추가할 때마다 [math(E)]의 핵을 계산한다. 즉, 다음을 만족하는 벡터 [math(\vec{a})]를 찾는다.

• [math(E\vec{a}\equiv\mathbf{0}\ (\mathrm{mod}\ 2))], [math(a_k=1)]

가우스 소거법을 이용해 핵을 효율적으로 관리할 수 있다. 만약 [math(\vec{a})]가 존재하지 않으면 [math(x_k)]를 하나 더 찾아야 한다.
마지막으로

• [math(\displaystyle x=\prod_{i=1}^kx^{a_i})]
• [math(\displaystyle y=\prod_{i=1}^kQ(x_i)^{a_i/2}=\prod_{j=0}^{\pi(B)}e_j^{\frac{1}{2}\sum_{i=1}^ka_ie_{i,\ j}})]

를 계산한다. [math(x^2\equiv y^2\ (\mathrm{mod}\ N))]이 성립함은 가정에서 드러나 있고, [math(x\not\equiv\pm y\ (\mathrm{mod}\ N))]이 성립하는지 확인하면 된다. 여기까지 왔다면, [math(\gcd(x-y,\ N))]을 계산했을 때 [math(N)]의 자명하지 않은 약수가 튀어나올 것이다. 안 그러면.. 수열에 원소를 추가하는 부분으로 다시 돌아가야 한다.

이차 체 알고리즘을 유사 코드로 나타내면 다음과 같다. 위의 설명과 달리 수열과 행렬이 0-based이기 때문에 변수의 index가 다른 점이 있다.

1. 홀수 [math(N)]을 입력받는다.
2. 소수 기저의 상한 [math(B)]를 정하고, 에라토스테네스의 체를 이용해 [math(B)] 이하의 소수를 모두 구하여 수열 [math(\left\{p_j\right\})]에 저장한다.
3. 빈 수열 [math(\left\{x_i\right\})], [math(\left\{A_i\right\})]와 크기가 [math(0\times(\pi(B)+1))]인 행렬 [math(E)], [math(V)], 크기가 [math(\pi(B)+1)]인 수열 [math(\left\{e_i\right\})]를 선언한다.[2]
4. [math(n=0)]
4. [math(x=\left\lfloor\sqrt{N}\right\rfloor)]
5. [math(Q(x)=0)]일 경우 [math(N)]은 제곱수이다. [math(x)]를 출력하고 종료한다.
6. [math(\displaystyle Q(x)=-C\prod_{j=0}^{pi(B)-1}p_j^{e_j})]와 같이 소인수분해한다. ([math(C)]는 [math(B)]보다 큰 소수의 곱으로 나타나는 자연수)
7. [math(C>1)]일 경우 11로
8. [math(e_{\pi(B)}=1)]
9. [math(x_0=x)], [math(A_0=1)], [math(E_0=e)], [math(V_0=e\mod2)]
10. [math(n\ \leftarrow\ n+1)]
11. [math(x^+=x)], [math(x^-=x)]
12. [math(x^+\ \leftarrow\ x^++1)]
13. [math(\displaystyle Q(x^+)=C\prod_{j=0}^{\pi(B)-1}p_j^{e_j})]와 같이 소인수분해한다.
14. [math(C>1)]일 경우 29로
15. [math(e_{\pi(B)}=0)], [math(a=2^n)], [math(\displaystyle v=\sum_{j=0}^{\pi(B)}2^j(e_j\mod2))]
16. [math(0\leq i<n)]인 정수 [math(i)]에 대해 다음을 반복한다.
17. [math(v=0)]일 경우, 25로
18. [math(V_n\ \leftarrow\ v)], [math(A_n\ \leftarrow\ a)]
19. [math(k=n)]
20. [math(k>0)]과 [math(V_k>V_{k-1})]이 모두 성립하는 동안, 다음을 반복한다.
21. [math(0\leq i<k)]인 정수 [math(i)]에 대해 다음을 반복한다.
22. [math(x_n\ \leftarrow\ x^+)], [math(E_n\ \leftarrow\ e)]
23. [math(n\ \leftarrow\ n+1)]
24. 29로
25. [math(y_1=x^+)]
26. [math(0\leq i<n)]인 정수 [math(i)]에 대해 다음을 반복한다.
27. [math(\displaystyle y_2=\prod_{j=0}^{\pi(B)-1}p_j^{\frac{e_j}{2}}\mod N)]
28. [math(y_1\neq y_2)]이고 [math(y_1+y_2\neq N)]일 경우, [math(\gcd(y_1-y_2,\ N))]을 출력하고 알고리즘을 종료한다.
29. [math(x^-\ \leftarrow\ x^--1)]
30. [math(\displaystyle Q(x^-)=C\prod_{j=0}^{\pi(B)-1}p_j^{e_j})]와 같이 소인수분해한다.
31. [math(C>1)]일 경우 46으로
32. [math(e_{\pi(B)}=0)], [math(a=2^n)], [math(\displaystyle v=\sum_{j=0}^{\pi(B)}2^j(e_j\mod2))]
33. [math(0\leq i<n)]인 정수 [math(i)]에 대해 다음을 반복한다.
34. [math(v=0)]일 경우, 42로
35. [math(V_n\ \leftarrow\ v)], [math(A_n\ \leftarrow\ a)]
36. [math(k=n)]
37. [math(k>0)]과 [math(V_k>V_{k-1})]이 모두 성립하는 동안, 다음을 반복한다.
38. [math(0\leq i<k)]인 정수 [math(i)]에 대해 다음을 반복한다.
39. [math(x_n\ \leftarrow\ x^+)], [math(E_n\ \leftarrow\ e)]
40. [math(n\ \leftarrow\ n+1)]
41. 46으로
42. [math(y_1=x^-)]
43. [math(0\leq i<n)]인 정수 [math(i)]에 대해 다음을 반복한다.
44. [math(\displaystyle y_2=\prod_{j=0}^{\pi(B)-1}p_j^{\frac{e_j}{2}}\mod N)]
45. [math(y_1\neq y_2)]이고 [math(y_1+y_2\neq N)]일 경우, [math(\gcd(y_1-y_2,\ N))]을 출력하고 알고리즘을 종료한다.
46. 12로

2.5.1. 최적화[편집]

---

여기서부터는 알고리즘 자체의 난이도가 크게 상승하는 대신, 수많은 최적화 기법을 통해 실행 속도를 크게 향상시킬 수 있다. 달리 말하면, 어떠한 최적화도 하지 않을 경우 폴라드 로 알고리즘보다도 느린 게 이차 체 알고리즘이다. 최적화에는 다음과 같은 요소가 있다.

• 오일러 조건을 이용한 소수 기저의 선별

• 토넬리-샹크스 알고리즘을 이용한 [math(x^2-N)] 소인수분해의 배치 처리

• 방정식 [math(x^2-N\equiv0\ (\mathrm{mod}\ p^n))]의 해의 점화식을 이용한 나눗셈 연산 횟수 줄이기

• [math(p)]가 홀수 소수이고 [math(x_n)]이 [math(x_n^2-N\equiv0\ (\mathrm{mod}\ p^n))]를 만족할 때, [math(x_{n+1}=x_n+(x_n^2+N)(-2x_1)^{-1})]이다. 이때 [math((-2x_1)^{-1})]은 [math(-2x_1)]의 mod [math(p)]에 대한 역원이다.
• [math(p=2)]일 때는 상황이 다르다. 먼저 [math(x_n^2-N\equiv0\ (\mathrm{mod}\ 2))]의 해는 모든 짝수이다. [math(N\equiv3\ (\mathrm{mod}\ 4))]이면 모든 짝수 [math(x)]에 대해 [math(x_n^2-N)]은 4의 배수가 아니고, [math(N\equiv5\ (\mathrm{mod}\ 8))]이면 모든 짝수 [math(x)]에 대해 [math(x_n^2-N)]은 8의 배수가 아닌 4의 배수이며, [math(N\equiv1\ (\mathrm{mod}\ 8))]이면 모든 짝수 [math(x)]에 대해 [math(x_n^2-N)]은 8의 배수이다. 16의 배수부터는 다음의 과정을 거친다.
  • [math(N\equiv1\ (\mathrm{mod}\ 16))]이면 [math(x_4=1)] 또는 [math(x_4=7)]이며,

• [math(n\ge4)]일 때, [math(x_{n+1}=x_n+(x_n^2-N)/2)]으로 계산한다.

• 희소 행렬을 이용한 [math(E)]와 [math(V)]의 메모리 최적화

2.6. 렌스트라 타원곡선 알고리즘[편집]

---

2.7. 다중 이차 체[편집]

---

2.8. 수체 체[편집]

---

2.9. 쇼어 알고리즘[편집]

---

양자컴퓨터에서 동작하는 알고리즘으로, 큐비트가 충분하다면 다항 시간에 소인수분해가 가능한 알고리즘. 자세한 내용은 문서 참조.