====== 소인수분해 (Prime Factorization) ====== * Integer Factorization 이라고도 불린다 (이게 더 흔히 쓰이는 표현같기도..) ===== n을 소인수분해하기 ===== * 참고 - [[https://cp-algorithms.com/algebra/factorization.html]] [[https://aruz.tistory.com/entry/factorization-1]] * 소인수분해는 암호학이나 여러 분야에서 매우 중요하게 사용되는 알고리즘이고, 그러다보니 다양한 알고리즘이 있는데.. PS 범위에서는 trial division과 Pollard's rho 두가지만 알면 충분하다 * 실제로는 더 효율적인 알고리즘들이 있지만.. 2^64 정도의 범위에서는 Pollard's rho 방법으로도 충분히 잘 돌아가는것 같다. * [[https://cp-algorithms.com/algebra/factorization.html#fermats-factorization-method|페르마의 방법]] 이나 [[https://cp-algorithms.com/algebra/factorization.html#pollards-p-1-method|Pollard's p-1]] 등은 Pollard's rho보다 장점이 없기 때문에 굳이 알 필요가 없지만, 소인수분해에 대한 이해를 위해서라면 읽어봐도 좋다 * [[https://koosaga.com/239]]은 더 효율적일수도 있지만, 복잡해서 PS에서는 못쓸것 같은 알고리즘이다. 한번 봐둬도 도움은 될듯. * 소인수분해 결과는 p_1^k_1 * p_2^k_2 * ... * p^n^k_n 의 형태로 나온다. 코드를 구현할때는 리턴값으로 p_i를 키로 k_i를 밸류로 갖는 dict를 리턴하도록 만들었다. n이 1일때는 빈 딕셔너리를 리턴한다 ==== Trial division ==== * n이 적당히 작은 경우에 효율적이다. * 그냥 i=2부터 시작해서 i*i def pollard_rho(n, x0, c): x, y = x0, (x0 * x0 + c) % n while (g := math.gcd(x - y, n)) == 1: x, y = (x * x + c) % n, (y * y + c) % n y = (y * y + c) % n return g * cpp 계열의 gcd 함수는 두 인자가 모두 양수여야 해서, abs(x-y)처럼 써줘야 하지만, python의 math.gcd는 abs를 안붙여도 제대로 된 값을 준다. * cpp에서는 x*x를 계산할때 오버플로우되지 않게 하기 위해서 추가 구현이 필요하지만, 파이썬에서는 신경쓸 필요가 없다. * x0와 c는 랜덤으로 생성된 값을 주거나, 그냥 고정된 값을 주는 경우도 있다. 고정된 값을 쓴다면 c=1, x0=2 를 쓰는것이 일반적이다 * 하지만 이 함수를 그대로 쓰기에는 문제가 있다. * 우선, 운나쁘게도 찾아낸 인수 g가 n과 같은 경우가 있다. 이 경우에는 f와 x0를 바꿔서 다시 돌리는 것이 정석이다. f(x)=x*x+c 의 형태는 유지한 채로, c값과 x0를 바꿔서 돌리면 된다. 랜덤을 다시 돌려서 찾거나, 이전값에서 1을 증가시키던가 해서 처리한다. * 이때 주의할점은, x0만 바꿔서는 돌려서는 안된다. 대표적인 예가 124376107291 을 인수분해 하는 경우인데, f=x*x+1 을 사용해서는 초기값을 아무리 바꿔도 인수를 찾아내지 못하고, 무한루프로 빠지게 된다. 차라리, x0는 안바꾸고 c만 바꾸는 것은 동작한다. * [[https://lpha-z.hatenablog.com/entry/2023/01/15/231500]] 를 참고하라 * [[https://github.com/cheran-senthil/PyRival/blob/master/pyrival/algebra/factors.py|pyrival]] 에서도 f를 고정시키고 x0만 바꾸는 방식의 잘못된 구현을 사용하고 있다. * 다음으로는 n이 소수일 경우이다. 이때는 인수가 진짜로 없으므로 f와 x0를 아무리 바꿔서 다시 돌려도 당연히 인수를 찾지 못한다. (g가 n으로 나온다). 따라서 폴라드 로를 돌리기 전에 미리 n이 소수인지 [[ps:소수 판별]]을 하고 소수가 아닌 경우에만 폴라드로를 돌리도록 한다. * 그래서, 폴라드로 알고리즘의 난이도는 폴라드로 자체보다도 오히려 소수판별을 위해 쓰는 밀러라빈의 난이도가 더 높다는 의견도 있다. * 지금처럼 밀러라빈+폴라드로를 세트로 묶어서 사용하는 것이 일반화되기 이전에는, 소수판별을 미리 하지 않고 그냥 돌려서 g==n일 경우 f와 x0를 바꿔서 몇번 더 돌려보고 계속 g==n이 나올경우에는 소수라고 간주하는 방식으로 돌리기도 했다고 한다. * 이부분을 추가하면 다음과 같아진다 * def find_factor(n): if is_prime(n): return n x0 = 2 for c in itertools.count(1): g = pollard_rho(n, x0, c) if g != n: return g * 소인수분해를 하기 위해서는 폴라드로 알고리즘으로 한개의 인수를 찾는것을 반복적으로 수행해야 한다 * n의 인수 f를 찾고 나면, f와 n/f에 대해서 다시 폴라드로를 돌려야 한다. 가장 간단한 것은 재귀로 구현하는 것이다. 아래는 pyrival의 구현체이다 * def prime_factors(n): """returns a Counter of the prime factorization of n""" if n <= 1: return Counter() f = find_factor(n) return Counter([n]) if f == n else prime_factors(f) + prime_factors(n // f) * pollard_rho(n, x0, c) 함수를 brent의 방법으로 바꿔서 속도를 더 빠르게 할수있다. * [[https://cp-algorithms.com/algebra/factorization.html#brents-algorithm]] 의 코드를 거의 그대로 python으로 옮기면 아래처럼 된다. * def brent_pollard(n,x0,c): m, l, x, g, q = 128, 1, x0, 1, 1 while g == 1: y = x for _ in range(l - 1): x = (x * x + c) % n k = 0 while k < l and g == 1: xs = x for _ in range(min(m, l - k)): x = (x * x + c) % n q = q * abs(y - x) % n g = math.gcd(q, n) k += m l += l if g != n: return g g = 1 while g == 1: xs = (xs * xs + 1) % n g = math.gcd(xs - y, n) return g * m 값에 대해서는 여기에서는 128을 줬는데, 1000을 사용하는 코드도 있다. * [[https://qiita.com/t_fuki/items/7cd50de54d3c5d063b4a]] 에서는 좀더 python에 맞춰 정리된 brent 알고리즘 코드가 있다. m을 n^(1/8)로 잡고 있다. * [[https://judge.yosupo.jp/submission/124476]]에서는 이 코드를 다시 깔끔하게 최적화한 코드를 볼수 있다. === 코드 === * [[:ps:teflib:numtheory#prime_factorization|teflib.numtheory.prime_factorization]] * 폴라드로 + brent 구현부분은 [[https://judge.yosupo.jp/submission/124476]]를 많이 참고했다 * 폴라드로를 반복적으로 돌리는 부분은 재귀 대신 반복으로 구현했다 === 관련 문제 === * ... ===== 여러개의 수를 소인수분해하기 ===== * x_1, x_2, ..., x_k 를 모두 소인수분해하기. 모든 수는 n이하라고 생각하자. * 소인수분해해야 할 수의 갯수 k에 따라서 다음중 적당한 방법들을 선택하면 된다 (n은 trial division으로 소인수분해가 가능한 범위라고 새생각하자) * k가 작으면, 그냥 위의 trial division 방법을 모든 x_i에 대해서 적용하면 된다. O(k*n^1/2) * k가 적당히 크면, 먼저 sqrt(n)이하의 소수 목록을 모두 계산해놓고, 그것에 대해서만 trial division을 사용하는 게 좀더 빠를수 있다 * 에라토스테네스의 체를 이용해서 소수 목록을 구하는 데에 O(sqrt(n)*loglogn) * 목록의 소수들만 이용해서 소인수분해를 하는데 O(sqrt(n)/ln(n)) * k개를 처리하는 총 시간복잡도는 대략 O(sqrt(n)*(loglogn + k/logn)). * 관련 문제 - [[ps:problems:boj:2312]] * k가 크면 (n에 가까울 정도로), 아래의 방법을 사용하자 O(n + klogn) ===== [1,n]까지의 모든 수를 소인수분해하기 ===== * 이 경우는, 모든 수에 대해서 최소인수를 미리 전처리해놓는 방식으로 구현한다. * 모든 수에 대해서 최소인수를 구하는 것은 보통 [[ps:소수 목록 구하기#linear sieve]]를 이용하는 것으로 설명하는 자료가 많은데, [[ps:소수 목록 구하기#에라토스테네스의 체]]를 써서도 가능하긴 하다. 하지만 [[ps:소수 목록 구하기#linear sieve]]를 이용하는 것이 조금 더 빠르다. 시간복잡도는 O(n). * 모든수의 최소인수를 구하고 나면, 임의의 수 x의 소인수분해를 O(Ω(x))에 처리할 수 있다. Ω(x)([[wp>Prime omega function]])은 x 소인수의 개수를 의미하는데, 최악의 경우 O(logx), 평균적으로는 O(loglogx)이다. * 결과적으로 1~n까지의 모든 수를 소인수분해하는 것도 O(nloglogn) 에 가능하다 * 코드 * [[:ps:teflib:numtheory#PrimeFactorizationCalculator|teflib.numtheory.PrimeFactorizationCalculator]] * 관련 문제 * [[ps:problems:boj:16563]] * [[ps:problems:boj:3142]]