[소인수분해] 폴라드 로 알고리즘 (Pollard's rho)

밀러-라빈 알고리즘(https://uneducatedjungler.tistory.com/121)에서 이어지는 글

 

폴라드 로 알고리즘은 입력된 정수 n의 랜덤한 인수를 뽑아내는 알고리즘이다.

정리에 앞서 거칠게 폴라드 로를 요약하자면 아래와 같다.

"n과 d의 공약수를 발견할 때까지 d에 유사난수를 대입한다"

 

얼핏 보면 n에 1부터 차례대로 나눠보는 것과 다름없는 굉장히 무식해보이는 방법같아 보이지만, 실제로는 의외로 빠르다. 그 이유는 생일 역설을 보면 알 수 있다.

 

생일 역설 (https://ko.wikipedia.org/wiki/%EC%83%9D%EC%9D%BC_%EB%AC%B8%EC%A0%9C)

 

생일은 366일이라는 한정된 range를 갖고 있다. 따라서 367명 이상이 있으면 반드시 생일이 같은 사람이 한 쌍은 존재하게 된다. 하지만 확률적으로 접근하면, 굳이 367명까지 가지 않더라도 60명 정도만 있어도 99.4%의 확률로 생일이 같은 한 쌍이 존재한다.

즉, n<367명의 사람이 있을 때 생일이 같은 한 쌍이 있을 확률은 우리가 인식하는 것보다 굉장히 높다는 것을 알 수 있다.

이를 다르게 표현해보면 아래와 같다.

모르는 사람의 생일을 맞추려 할 때, 랜덤한 생일을 아무렇게나 찔러봐도 한 60번정도면 맞출 확률이 99.4%나 된다.

 

여기서 60이라는 try의 수와 99.4%라는 확률은 어떻게 나온 걸까? 그것은 366이라는 range의 값으로부터 유래한다. n번의 try 안에 찾을 확률 p는 (364/365)*(363/365)* ... *((365-n)/365)인데, 이 값이 range가 커질수록 표준편차 root(range)의 정규분포에 수렴하기 때문이다.

(366의 제곱근은 약 19.1이고, 60은 그 약 3배의 값이다. 그리고 정규분포에서 3시그마 이내에 값이 존재할 확률은 99.5%이다.)

 

일반화하면, 유한한 range를 가진 값에서 중복값을 발견할 확률은 sigma = root(range)의 정규분포를 따른다. 다시 말해, 대충 2*root(range)번만 비교해봐도 95%의 확률로 발견할 수 있고, 3*root(range)번 비교해보면 99.5%의 확률로 발견이 가능하다는 것이다.

 

폴라드 로 알고리즘은 바로 이점을 이용하여 평균 O(root(n))의 시간복잡도로 n의 인수 하나를 발견한다.

위의 볼드한 문장과 같은 형식으로 써보자면, 아래와 같은 의미이다.

모르는 수의 인수를 맞추려 할 때, 랜덤한 수를 아무렇게나 찔러봐도 한 3root(n)번정도면 맞출 확률이 99.4%나 된다.

 

import sys
input = sys.stdin.readline
sys.setrecursionlimit(10**8)
from math import gcd
from random import randrange
 
def MillerRabin(n) : 
    // 구현 생략
    // n이 소수일시 return True, 아닐시 return False
 
def f(x, c, n) : 
    return (x**2 + c) % n
 
def PollardRho(n) : 
    if n ==  1 : 
        return 1
    if n % 2 == 0 : 
        return 2
    if MillerRabin(n) : 
        return n
    x = randrange(2, n)
    c = randrange(1, n)
    y = x
    g = 1
    while g == 1 : 
        x = f(x, c, n)
        y = f(f(y, c, n), c, n)
        g = gcd(x - y, n)
        if g == n : 
            return PollardRho(n)
    if MillerRabin(g) : 
        return g
    else : 
        return PollardRho(g)
 
n = int(input())
factor = []
while n != 1 : 
    p = PollardRho(n)
    n = n // p
    if not p in factor : 
        factor.append(p)

일정한 규칙으로 생성된 유사난수가 있을 때, 이를 mod n한 값은 0~n-1중 하나이다. (= range가 n이다.)

즉 규칙적으로 유사난수 생성 시 최악의 경우에도 n번 이내에 사이클이 발생하여 최초의 값으로 돌아오게 된다. 이는 최악의 경우에도 while문이 무한히 루프하지 않고 언젠가 함수가 종료됨을 보장한다. (유사난수를 쓰는 이유) 일반적으로는 x^2+1 등 효율성이 잘 알려진 유사난수 패턴을 이용한다. (line 11 ~ 12, 21, 22)

내가 구현한 위 라이브러리에서는 x^2+1로 추출 불가능한 경우에 대해 다른 패턴을 대입하기 위해 상수항 c에도 randrange를 썼다.

 

 

이제 이 유사난수와 n의 공약수가 1이나 n이 아닐 때까지 while문을 돌려가며 호제법을 쓰면, root(n) 단위의 시간에 인수 하나를 발견할 수 있는 것이다. (line 21 ~ 30)

여기까지가 랜덤한 인수를 발견하는 폴라드 로 알고리즘의 기초이고, 이를 빠른 속도의 소수판별법(일반적으로 밀러-라빈)과 같이 응용하면 n을 완전히 소인수분해하는 것도 가능하다.

 

폴라드 로로 발견한 인수 g에 대해, 앞서 정리한 밀러 라빈 소수판별법을 이용하여 빠르게 g가 소수인지 아닌지 판별하여, 소수라면 return하고, 소수가 아니라면 g에 대해 다시 폴라드 로 알고리즘을 실행한다. (line 31 ~ 34)

이러한 재귀반복을 통해 n의 한 소인수 p를 발견할 수 있다.

 

이제 n = n//p에 대해서도 위의 과정을 반복하며, n이 1이 될때까지 소인수분해를 계속하면 n을 완전히 소인수분해할 수 있다. (line 36 ~ 42)

 

 

폴라드 로를 PS에 응용해서 풀 수 있는 문제의 예는 다음과 같다.

4149 큰 수 소인수분해 (https://www.acmicpc.net/problem/4149)

17633 제곱수의 합 (More Huge) (https://www.acmicpc.net/problem/17633)

 

위 라이브러리는 정말 소인수만을 구하기 위한 코드라 약간의 비효율이 있다.

예를 들어 line 39에서 추출한 p에 대해 나머지가 0인 동안 계속 나누는 것으로 폴라드로를 다시한번 돌렸는데 같은 소수가 나와서 시간을 낭비하는 경우를 줄일 수 있다.

그런데 이때 while문이 도는 동안 order를 +1하고, factor에 [p, order]를 담는다면 소인수를 그 차수와 함께 구할 수 있을 것이다. 이 [소인수, 차수] 쌍은 오일러 피 함수의 값을 구하는 데에도 사용할 수 있다.