Project Euler 3: 나눌 수 있는 만큼 나눠서 가장 큰 소인수 찾기

스포일러 주의 · 이 글에는 풀이 과정과 정답이 그대로 나옵니다. 직접 먼저 풀어 보고 싶다면 Project Euler 3번 문제부터 보세요.

Project Euler 시리즈의 세 번째 문제입니다. 소인수분해는 이 시리즈에서 계속 다시 나오는 도구인데, 그 가장 기본적인 형태인 시행 나눗셈(trial division)을 이번에 정확히 다져 둡니다. 코드는 열 줄 남짓이지만 “반복이 끝났을 때 남은 값이 왜 소수인가”, “반복은 최악의 경우 몇 번 도는가”를 제대로 따져 보면 생각보다 배울 것이 많습니다.

스포일러 주의: 이 글에는 풀이 전체와 정답이 들어 있습니다. 직접 먼저 풀어 보고 싶다면 원래 문제부터 보세요. 이 시리즈의 글에는 광고를 싣지 않습니다. 문제 본문의 라이선스(CC BY-NC-SA 4.0)가 비영리 조건을 두기 때문입니다.

문제

어떤 정수를 나누어떨어지게 하는 소수를 그 수의 소인수라고 합니다. 문제에 나온 예로 13195의 소인수는 5, 7, 13, 29이고, 그중 가장 큰 것은 29입니다. 같은 방식으로 600851475143의 가장 큰 소인수를 구하는 것이 문제입니다.

입력은 12자리 정수 하나이고, 구할 것은 소수 하나입니다.

필요한 배경 지식

소수와 합성수. 2 이상의 정수 중 약수가 1과 자기 자신뿐인 수가 소수, 그 밖의 2 이상의 정수가 합성수입니다. 1은 소수도 합성수도 아니고 소인수가 하나도 없습니다.

소인수분해와 유클리드의 보조정리. 2 이상의 모든 정수는 소수의 곱으로 나타낼 수 있고, 곱하는 순서를 빼면 그 방법은 하나뿐입니다(산술의 기본 정리). 이 글에서 직접 쓰는 것은 그 바탕이 되는 유클리드의 보조정리입니다. 소수 p가 곱 ab를 나누면 p는 a나 b 중 적어도 하나를 나눕니다. 그래서 n = a × b로 쪼개 두면, n의 소인수는 a의 소인수이거나 b의 소인수입니다.

합성수에는 √m 이하의 소인수가 있습니다. m이 합성수이면 m = a × b (1 < a ≤ b < m)로 쓸 수 있습니다. a ≤ b이므로 a × a ≤ a × b = m, 즉 a ≤ √m입니다. a의 소인수 하나를 p라고 하면 p ≤ a ≤ √m이고 p는 m도 나눕니다. 대우로 말하면, 2 이상의 m이 √m 이하의 어떤 소수로도 나누어떨어지지 않으면 m은 소수입니다. 이 사실이 이번 풀이의 핵심입니다.

정수 나눗셈과 내림. 양의 정수 f, r에 대해 r / f를 정수 나눗셈(소수점 아래 버림, ⌊r / f⌋)으로 계산할 때, f × f ≤ r과 f ≤ ⌊r / f⌋는 같은 조건입니다. f × f ≤ r이면 f ≤ r / f이고, f는 정수이므로 r / f를 내린 값 이하이기도 합니다. 거꾸로 f ≤ ⌊r / f⌋ ≤ r / f이면 양변에 f를 곱해 f × f ≤ r입니다. 이 동치는 C++ 구현에서 곱셈 넘침을 피하는 데 씁니다.

가장 단순한 방법

정의를 그대로 옮기면 2부터 n까지 모든 정수 d를 보면서 “d가 n을 나누는가”, “d가 소수인가”를 확인하고, 둘 다 맞는 d 중 가장 큰 것을 고르면 됩니다.

def is_prime(m: int) -> bool:
    return m >= 2 and all(m % d != 0 for d in range(2, m))

n = 13195
best = 0
for d in range(2, n + 1):
    if n % d == 0 and is_prime(d):
        best = d

작은 예에서는 문제없이 29가 나오지만, n이 600851475143이면 바깥 반복만 6천억 번 가까이 돌아야 합니다. 바깥 반복만 해도 O(n)이고, 약수마다 소수 판정을 따로 하면 더 듭니다. n의 자릿수가 하나 늘 때마다 일이 열 배씩 늘어나므로, 이 방법은 n이 커지면 쓸 수 없습니다.

위쪽에서 아래로 내려오며 처음 만나는 소수 약수에서 멈추는 식으로 바꿔도, 가장 큰 소인수가 n보다 훨씬 작으면 그사이의 수를 거의 모두 확인해야 하므로 근본적으로 나아지지 않습니다.

풀이 유도

1단계: 작은 소인수를 나눠 없애면 문제가 작아집니다. 13195를 작은 수부터 나눠 봅니다. 2로는 안 나누어떨어지고(홀수), 3도 안 됩니다(자릿수 합 19). 4는 2가 안 되니 당연히 안 됩니다. 5로 나누면

13195 = 5 × 2639

입니다. 이제 13195의 나머지 소인수는 모두 2639의 소인수입니다(유클리드의 보조정리). 원래 수 대신 더 작은 2639만 보면 됩니다.

2단계: 같은 일을 이어 갑니다. 2639는 5로 더 나누어떨어지지 않으므로 6으로 넘어가고, 6도 안 되니 7로 나눕니다.

2639 = 7 × 377

377은 7로 더 나누어떨어지지 않습니다. 8, 9, 10, 11, 12는 모두 안 되고, 13에서

377 = 13 × 29

가 됩니다. 남은 값은 29입니다.

3단계: 언제 멈출 수 있는지 봅니다. 다음 후보는 14인데 14 × 14 = 196 > 29입니다. 지금까지 2부터 13까지의 수로 나눌 수 있는 만큼 다 나눴으므로 29는 13 이하의 어떤 소수로도 나누어떨어지지 않습니다. 그런데 √29 ≈ 5.4이므로 √29 이하의 소수는 이미 모두 확인했습니다. 배경 지식의 사실에 따라 29는 소수입니다. 29를 더 시험해 볼 필요가 없습니다.

정리하면 13195 = 5 × 7 × 13 × 29이고, 나눠 없앤 5, 7, 13보다 마지막에 남은 29가 크므로 답은 29입니다. 문제의 예와 같습니다.

4단계: 합성수 후보는 걱정하지 않아도 됩니다. 1단계에서 4나 6 같은 합성수도 후보로 시험했지만 나누어떨어지지 않았습니다. 우연이 아닙니다. 후보 f가 합성수이면 f의 소인수는 모두 f보다 작고, 그 소수들은 f 차례가 오기 전에 이미 남은 값에서 다 나눠 없앴습니다. 그러니 남은 값이 f로 나누어떨어질 수 없습니다. 결과적으로 실제로 나누어떨어지는 후보는 언제나 소수입니다. 소수만 골라 시험하려고 소수 목록을 따로 만들 필요가 없습니다.

최종 알고리즘

n < 2이면 소인수가 없으므로 0을 돌려줍니다(아래 “경계 사례”).
remaining = n, largest = 1, factor = 2로 시작합니다.
factor × factor ≤ remaining인 동안, remaining이 factor로 나누어떨어지는 한 계속 나누고 그때마다 largest = factor로 둡니다. 그다음 factor를 1 늘립니다.
반복이 끝났을 때 remaining > 1이면 largest = remaining입니다.
largest를 돌려줍니다.
정답을 출력하기 전에 작은 입력을 확인합니다. 1 → 0, 2 → 2, 9 → 3, 64 → 2, 97 → 97, 13195 → 29. 하나라도 다르면 0이 아닌 종료 코드로 끝냅니다.

C++ 구현은 3단계의 조건을 factor ≤ remaining / factor로 씁니다. 배경 지식에서 본 대로 같은 조건이지만 곱셈이 없어 넘치지 않습니다. 이유는 “정수 범위”에서 설명합니다.

정확성 증명

반복문이 후보 f로 조건을 검사하는 순간(f = 2, 3, …)마다 다음이 성립합니다(불변식).

(가) n = (지금까지 나눠 없앤 소수들의 곱) × remaining이고, 나눠 없앤 소수는 모두 f보다 작습니다.
(나) remaining은 f보다 작은 어떤 소수로도 나누어떨어지지 않습니다.
(다) 하나라도 나눠 없앴다면 largest는 그중 가장 큰 것이고, 아니면 1입니다.

처음(f = 2)에는 나눠 없앤 것이 없고 remaining = n, largest = 1이므로 성립합니다. 조건이 참이어서 f 차례를 처리할 때, 안쪽 반복은 remaining이 f로 나누어떨어지지 않을 때까지 f를 나눠 냅니다. f가 실제로 나누었다면 4단계에서 본 대로 f는 소수이고, 나눌 때마다 곱에 f를 하나 옮기므로 (가)가 유지됩니다. 안쪽 반복이 끝나면 remaining은 f로 나누어떨어지지 않으므로 f + 1에서 (나)가 성립합니다. f는 지금까지 나눠 없앤 어떤 소수보다 크므로 largest = f로 두면 (다)도 유지됩니다. remaining은 한 번 나눌 때마다 적어도 절반이 되므로 안쪽 반복은 반드시 끝납니다.

바깥 반복은 f가 1씩 늘고 remaining은 줄지 않으므로 언젠가 f × f > remaining이 되어 끝납니다. 끝나는 순간의 f를 F라고 하면 두 경우가 있습니다.

remaining = 1: (가)에서 n은 나눠 없앤 소수들의 곱입니다. 유클리드의 보조정리에 따라 n의 소인수는 이 소수들 중 하나이므로, (다)의 largest가 가장 큰 소인수입니다.
remaining > 1: (나)에서 remaining은 F보다 작은 소수로 나누어떨어지지 않고, F × F > remaining이므로 √remaining 이하의 소수는 모두 F보다 작습니다. 배경 지식의 사실에 따라 remaining은 소수입니다. 또 remaining의 가장 작은 소인수(=remaining 자신)는 F 이상이므로, F보다 작은 나눠 없앤 소수들보다 큽니다. n의 소인수는 나눠 없앤 소수이거나 remaining이므로 가장 큰 소인수는 remaining입니다.

두 경우 모두 반환값은 n의 가장 큰 소인수입니다.

복잡도

최악의 반복 횟수. 바깥 반복의 본문은 f × f ≤ remaining ≤ n인 f에 대해서만 돌므로 많아야 ⌊√n⌋ − 1번입니다. 안쪽 반복에서 실제로 나누는 횟수는 n을 소인수분해했을 때 소인수의 개수(중복 포함)와 같고, 소인수는 모두 2 이상이므로 log₂ n 이하입니다. 그 밖의 연산은 바깥 반복 한 번마다 나머지 검사 한 번과 조건 검사 한 번입니다. 고정 폭 정수에서 이 연산들은 상수 시간이므로 시간 복잡도는 O(√n) 이고, 저장하는 값이 세 개뿐이라 추가 공간은 O(1) 입니다.

이 상한은 실제로 도달합니다. n이 소수이면 remaining이 한 번도 줄지 않아 f가 √n 근처까지 올라가야 합니다. “나누면서 remaining이 줄어드니 로그 시간”이라고 하면 틀립니다. 나누는 횟수는 로그로 제한되지만, 나누어떨어지지 않는 후보를 지나가는 횟수는 그렇지 않습니다.

입력에 따라 빨리 끝나는 경우. 반복은 remaining이 작아질수록 일찍 멈춥니다. 소인수를 중복을 포함해 크기순으로 늘어놓았을 때 가장 큰 것을 p₁, 두 번째를 p₂라고 하면, 반복은 대략 max(p₂, √p₁)까지만 돕니다. p₂까지는 나눠 내야 하고, 그 뒤에는 remaining이 p₁ 하나만 남아 f가 √p₁을 넘으면 멈추기 때문입니다. 작은 소인수가 많은 수일수록 √n보다 훨씬 일찍 끝납니다. 이 글의 입력에서 실제로 몇 번 도는지는 “정답과 검증”에서 셉니다.

정수 범위.

600851475143은 2^39 ≈ 5.5 × 10^11보다 크고 2^40 ≈ 1.1 × 10^12보다 작은 40비트 수입니다. 32비트 정수(최대 약 4.3 × 10^9)에는 들어가지 않으므로 64비트 정수가 필요합니다.
C++는 unsigned long long을 씁니다. C++ 표준은 이 타입이 최소 64비트임을 보장하고(C 표준의 ULLONG_MAX 최솟값을 따름), 부호 없는 연산은 2^N을 법으로 정의되어 있어 정의되지 않은 동작이 생기지 않습니다. 다만 넘치면 조용히 틀린 값이 됩니다. 리터럴에는 ULL 접미사를 붙여 처음부터 같은 타입으로 둡니다.
반복 안의 나눗셈과 나머지 연산은 값을 키우지 않으므로 넘칠 일이 없습니다. 위험한 것은 조건 factor * factor <= remaining의 곱셈입니다. remaining이 2^64에 가까운 소수라면 factor가 2^32까지 올라가고, 그때 factor × factor = 2^64는 0으로 돌아가 조건이 다시 참이 됩니다. 그러면 반복이 끝나야 할 때 끝나지 않습니다. 그래서 C++ 구현은 같은 조건을 factor <= remaining / factor로 씁니다. 나눗셈만 쓰므로 넘치지 않고, 배경 지식의 동치로 결과도 같습니다.
factor 자체는 조건이 참인 동안 √remaining < 2^32 이하이므로 ++factor도 넘치지 않습니다. 따라서 C++ 함수는 unsigned long long의 모든 값에서 정확합니다. 다만 2^64에 가까운 소수를 넣으면 반복이 약 2^32 ≈ 43억 번 돌아야 하므로 정확하더라도 느립니다.
Python의 int는 정밀도에 제한이 없어 넘칠 일이 없으므로 조건을 factor * factor <= remaining 그대로 썼습니다. 알고리즘은 같습니다.

경계 사례

n < 2: 0과 1은 소인수가 없어 “가장 큰 소인수”가 정의되지 않습니다. 두 구현 모두 이 경우 0을 돌려주기로 정했습니다. 0은 소수가 아니므로 정상적인 답과 헷갈리지 않습니다. Python에서는 음수도 같은 처리를 받고, C++는 부호 없는 타입이라 음수를 받지 않습니다. 작은 예 검사에 1 → 0이 들어 있습니다.
n이 소수: 반복 안에서 한 번도 나누지 못하고 끝나며, 남은 n 자신이 답입니다(2 → 2, 97 → 97). 2는 처음부터 2 × 2 > 2라서 반복 본문이 한 번도 돌지 않습니다.
n이 소수의 거듭제곱: 64 = 2^6은 2로 여섯 번 나누면 remaining이 1이 되고, 답은 마지막으로 나눈 2입니다. 마지막 단계의 remaining > 1 검사가 1을 답으로 착각하지 않게 막습니다.
n이 소수의 제곱: 9는 f = 3에서 3 × 3 ≤ 9로 조건이 참이어야 나눠 낼 수 있습니다. 조건을 <로 잘못 쓰면 9를 소수로 판단해 9를 돌려주게 되므로, 이 경우가 작은 예 검사에 들어 있습니다.

C++17 구현

#include <iostream>

// Largest prime factor of n. Values below 2 have no prime factor; return 0 for them.
// Divide out every factor from 2 upward while factor * factor <= remaining.
// The loop tests factor <= remaining / factor, the same condition without a
// multiplication that could overflow. Whatever is left above 1 at the end is
// prime and larger than every factor divided out before it.
unsigned long long largest_prime_factor(unsigned long long n) {
    if (n < 2) {
        return 0;
    }
    unsigned long long remaining = n;
    unsigned long long largest = 1;
    for (unsigned long long factor = 2; factor <= remaining / factor; ++factor) {
        while (remaining % factor == 0) {
            remaining /= factor;
            largest = factor;
        }
    }
    if (remaining > 1) {
        largest = remaining;
    }
    return largest;
}

int main() {
    if (largest_prime_factor(1) != 0 || largest_prime_factor(2) != 2 ||
        largest_prime_factor(9) != 3 || largest_prime_factor(64) != 2 ||
        largest_prime_factor(97) != 97 || largest_prime_factor(13195) != 29) {
        return 1;
    }
    std::cout << largest_prime_factor(600851475143ULL) << '\n';
    return 0;
}

작은 예의 확인은 assert가 아니라 if로 했습니다. assert는 NDEBUG가 정의된 빌드에서 사라지기 때문입니다. 확인이 하나라도 틀리면 return 1로 0이 아닌 종료 코드를 냅니다.

Python 3 구현

def largest_prime_factor(n: int) -> int:
    """Largest prime factor of n. Values below 2 have no prime factor; return 0 for them.

    Divide out every factor from 2 upward while factor * factor <= remaining.
    Whatever is left above 1 at the end is prime and larger than every factor
    divided out before it.
    """
    if n < 2:
        return 0
    remaining = n
    largest = 1
    factor = 2
    while factor * factor <= remaining:
        while remaining % factor == 0:
            remaining //= factor
            largest = factor
        factor += 1
    if remaining > 1:
        largest = remaining
    return largest


def main() -> None:
    if (
        largest_prime_factor(1) != 0
        or largest_prime_factor(2) != 2
        or largest_prime_factor(9) != 3
        or largest_prime_factor(64) != 2
        or largest_prime_factor(97) != 97
        or largest_prime_factor(13195) != 29
    ):
        raise SystemExit(1)
    print(largest_prime_factor(600851475143))


if __name__ == "__main__":
    main()

C++와 같은 알고리즘, 같은 작은 예입니다. 나눗셈은 정수 나눗셈 //=로 해서 값이 float로 바뀌지 않게 합니다. /=를 쓰면 결과가 부동소수점이 되어 큰 수에서 정밀도를 잃을 수 있습니다. 작은 예가 틀리면 raise SystemExit(1)로 종료 코드 1을 냅니다. assert는 python -O로 실행하면 제거되므로 쓰지 않았습니다.

정답과 검증

600851475143의 가장 큰 소인수는 6857입니다.

알고리즘이 지나가는 과정은 다음과 같습니다.

factor = 2 … 70      나누어떨어지지 않음
factor = 71          600851475143 = 71 × 8462696833
factor = 72 … 838    나누어떨어지지 않음
factor = 839         8462696833   = 839 × 10086647
factor = 840 … 1470  나누어떨어지지 않음
factor = 1471        10086647     = 1471 × 6857
factor = 1472        1472 × 1472 > 6857 → 멈춤, 남은 6857은 소수

즉 600851475143은 표에 나온 네 소수의 곱입니다. 바깥 반복의 본문은 factor = 2부터 1471까지 1470번 돌고, 실제 나눗셈은 세 번입니다. 처음 상한인 ⌊√600851475143⌋ = 775146보다 훨씬 일찍 끝나는 것은 “복잡도”에서 본 대로 두 번째로 큰 소인수가 1471로 작기 때문입니다.

검증은 세 겹으로 합니다.

작은 예. 두 프로그램 모두 답을 출력하기 전에 1, 2, 9, 64, 97, 13195의 결과가 0, 2, 3, 2, 97, 29인지 확인하고, 다르면 0이 아닌 종료 코드로 끝납니다. 13195 → 29는 “풀이 유도”에서 손으로도 따라가 보았습니다.
두 구현의 일치. 이 시리즈의 저장소에서는 검사 스크립트가 C++ 코드를 g++ -std=c++17 -O2 -Wall -Wextra -Werror로 컴파일해 실행하고(경고도 실패로 처리), Python 코드를 python3 -I -B로 실행합니다. 두 출력이 저장해 둔 answer.txt와 바이트 단위로 같아야 하고, 0이 아닌 종료 코드나 시간 초과(C++ 1초, Python 10초)는 실패입니다. 시리즈 규칙상 이 검사를 통과한 코드만 글에 싣습니다.
손 계산. 위 표의 네 소인수를 다시 곱하면 71 × 839 = 59569, 뒤의 두 수의 곱은 10086647이고, 59569 × 10086647 = 600851475143으로 원래 수가 나옵니다. 네 수가 소수인지도 배경 지식의 사실로 확인할 수 있습니다. 예를 들어 가장 큰 소인수는 제곱근이 약 82.8이고, 그 이하의 소수 2, 3, 5, …, 79 중 어느 것으로도 나누어떨어지지 않습니다.

더 생각해 볼 것

후보를 반으로 줄이기. 2를 먼저 다 나눠 낸 뒤에는 홀수 후보만 시험해도 됩니다. 정확성 증명의 어느 부분을 고쳐야 하는지, 반복 횟수가 얼마나 줄어드는지 따져 보세요. 2와 3을 먼저 처리하고 6k ± 1 꼴만 시험하면 어떨까요?
모든 소인수와 지수. 같은 반복에서 나눌 때마다 (소수, 횟수)를 기록하면 소인수분해 전체가 나옵니다. 거기서 약수의 개수와 합을 구하는 공식으로 이어지는데, 이 시리즈의 1–100번 안에서도 약수의 개수나 합을 다루는 문제(예: 12번, 21번)가 나옵니다.
많은 수를 한꺼번에. 1부터 N까지 모든 수의 가장 작은 소인수를 미리 체로 구해 두면, 각 수를 그 표를 따라 몇 번만 나눠서 분해할 수 있습니다. 에라토스테네스의 체를 어떻게 바꾸면 되는지 생각해 보세요. 소수를 많이 다루는 이 시리즈의 7번, 10번에서도 체가 쓸모 있습니다.
최악의 입력. 64비트 범위에서 이 코드가 가장 오래 걸리는 입력은 어떤 수일까요? 큰 소수와 비슷한 크기의 두 소수의 곱은 각각 언제 멈추는지 “복잡도”의 max(p₂, √p₁)로 비교해 보세요.
넘침 실험. C++ 조건을 factor * factor <= remaining로 바꾸고, 2^64보다 작은 가장 큰 소수 근처의 값을 넣으면 어떤 일이 생기는지 생각해 보세요(실행하면 오래 걸립니다).

출처와 라이선스

Project Euler, Problem 3 — 원래 문제. 이 글의 “문제” 절은 그 내용을 직접 요약한 것이고, 원문과 같은 CC BY-NC-SA 4.0으로 공유합니다. 문제의 저작자 표시는 Project Euler입니다.
Project Euler, Copyright — 문제 본문의 라이선스와 링크 방식. 비영리 조건에 따라 이 글에는 광고를 싣지 않습니다.
풀이, 증명, 코드는 직접 쓴 것입니다. Project Euler의 해설(overview), 포럼, 회원 전용 자료나 다른 풀이 사이트는 참고하지 않았습니다. Motitory는 Project Euler와 제휴나 후원 관계가 없고, Project Euler가 이 글을 검토하지도 않았습니다.
C++17 최종 작업 초안 N4659 — §6.9.1 기본 타입(부호 없는 정수의 2^N 법 연산), §21.3.6 <climits>(C 표준의 범위를 따름), §22.3 assert(NDEBUG)
Python 표준 라이브러리: Numeric Types — int의 무제한 정밀도, //(내림 나눗셈)과 /(실수 나눗셈)
Python 레퍼런스: The assert statement — 최적화 모드에서 assert가 제거됨