오일러 프로젝트 12
1부터 n까지의 자연수를 차례로 더하여 구해진 값을 삼각수라고 합니다. 예를 들어 7번째 삼각수는 1 + 2 + 3 + 4 + 5 + 6 + 7 = 28이 됩니다. 이런 식으로 삼각수를 구해 나가면 다음과 같습니다.
1, 3, 6, 10, 15, 21, 28, 36, 45, 55, ...이 삼각수들의 약수를 구해봅시다.
1: 1 3: 1, 3 6: 1, 2, 3, 6 10: 1, 2, 5, 10 15: 1, 3, 5, 15 21: 1, 3, 7, 21 28: 1, 2, 4, 7, 14, 28위에서 보듯이, 5개 이상의 약수를 갖는 첫번째 삼각수는 28입니다. 그러면 500개 이상의 약수를 갖는 가장 작은 삼각수는 얼마입니까
가장 작은 삼각수부터 약수의 개수를 계산하여, 처음으로 약수의 개수가 500을 넘는 값을 찾으면 된다.
삼각수
삼각수는 파스칼의 삼각형에 나오는 수이고, N번째 삼각수는 1부터 N까지의 합으로 계산된다.
\large A_{n} = \frac{n \times (n + 1)}{2}
i = 0, 1, 2, 3, 4, … 일 때 i를 위 공식에 대입해서 구해도 되지만, 1씩 증가하는 값을 누적으로 구해서 사용해도 된다. 만약 어떤 자연수의 약수의 개수를 세는 numbers_of_factors()
라는 함수를 우리가 갖고 있다면, 이 문제는 아래와 같은 코드로 풀어낼 수 있다.
%%time
a, b = 1, 2
while True:
if number_of_factors(a) > 500:
break
a, b = a + b, b + 1
print(a)
약수의 개수 세기
어떤 자연수 N이 자연수 p로 나눠진다면, 두 개의 약수가 있다는 것을 의미한다. 따라서 1부터 N의 거듭제곱까지의 범위에서 각각의 자연수로 나누고, 나눌 수 있다면 약수가 2개씩 있다고 생각하면 된다. 단, N이 완전제곱수인 경우에는 제곱근인 약수는 두 개로 세면 안된다. 이 부분만 주의해서 약수의 개수를 세는 함수를 아래와 같이 작성할 수 있다.
def numbers_of_factors(n):
l = int(n**0.5)
s = 1 if l * l == n else 2
a = 2
while a <= l:
if n % a == 0:
s += 2
a += 1
return s
최적화
전체 코드는 앞에서 언급했으니, 실행하면 답을 구할 수 있다. 단, 실제 정답은 꽤 큰 수이기 때문에 삼각수를 순서대로 만드는 방법을 사용하지 않고 1씩 증가시키면서 삼각수인지 검사하려 든다면, 왠만한 전기세로는 답을 볼 수가 없을 것이다. 또한 실제로 이 방법을 사용해도 시간은 꽤 오래 걸리는 편이다.
실제로 적당히 빠른 시간 (3초 이내?) 내에 답을 구하기 위해서는 어느 정도 최적화를 해줘야 하는 문제가 거의 처음으로 나왔다고 볼 수 있겠다. 약수의 개수를 세는 연산이, 비록 제곱근까지만 검사해서 전체 루프 수를 줄이기는 하지만, N값이 크면 클수록 분명 더 많은 루프를 돌면서 나눗셈을 해봐야 하는 것이다. 그래서 시간을 더 단축하고 싶다면, 약수의 개수를 더 빠르게 구하는 방법을 사용해야 한다.
소인수 분해하기
소인수분해는 2부터 1씩 증가하면서 n을 나눌 수 있을만큼 나누고, 나눈 수와 나눈 횟수를 기록해둔다. 이걸 한계값 (제곱근인데, N을 나눠서 더 작은 N’ 가 되었다면 N’의 제곱근까지만 검사하면 된다)까지 반복하고 나눌 수 있는 인수와 그 개수를 리턴해주면 된다.
N이 아주 크다면 이 마저도 ‘나눌 수 있는 다음 후보’를 똑똑하게 찾도록 튜닝을 해야하는데, 일단 최근에는 PC 성능도 좋아졌고, 파이썬 자체의 성능도 좋아졌기 때문에 이 문제에서는 이 정도 함수면 그래도 1초 이내에서 답을 구할 수 있다.
소인수분해를 하면 각각의 소인수와 그 지수들을 얻게 된다. 각 약수는 각각의 소인수를 인수로 포함하거나 하지 않을 수 있기 때문에, 이러한 경우를 모두 세어보면 각 인수의 지수에 + 1씩을 더한 후 모두 곱하면 된다. 따라서 소인수분해 함수를 사용하여 약수의 개수를 구하는 함수는 아래와 같이 작성할 수 있다.
from functools import reduce
def factorize(n):
res = dict()
e = 0
while n % 2 == 0:
n, e = n // 2, e + 1
if e > 0:
res[2] = e
k = 3
while k < int(n ** .5 + 1.5):
e = 0
while n % k == 0:
n, e = n // k, e + 1
if e > 0:
res[k] = e
k += 2
if n > 1:
res[n] = 1
return res
def numbers_of_factor(n):
if n == 1:
return 1
return reduce(lambda x, y: x * y, (x + 1 for x in factorize(n).values()))
%%time
a, b = 1, 2
while numbers_of_div(a) < 500:
a, b = a + b, b + 1
print(a)
# 76576500
# CPU times: totla: 750 ms