[Project Euler] 44. Pentagon numbers

1. 개요

    문제는 이곳에서 확인할 수 있다. pentagonal number \(p_n\)은 \(n(3n-1)/2\)로 정의되며 처음 10개는 다음과 같다.

1, 5, 12, 22, 35, 51, 70, 92, 117, 145, ...

    합과 차가 pentagonal number인 어떤 두 pentagonal number가 존재할때, 가장 작은 차 \(D=|P_k-P_j|\)를 구하는 것이 문제의 목표이다.

2. 구현

  합은 \(P_l=P_k+P_j\)로 차는 \(P_i=P_k-P_j\)로 정의하자. 또한 \(i<j<k<l\)이다. 가장 먼저 생각한 방법은 \(k\)를 순차적으로 늘려가면서 \(P_l\)과 \(P_i\)를 찾는 것이었다. 가장 작은 \(P_l\)에 대응하는 \(P_i\)가 최소 차 즉, \(D\)일 것이라고 생각했다. 또한 어떤 수 \(x\)가 pentagonal number인지 확인하는 방법으로 근의 공식을 사용하였다. 구현은 다음과 같다.

# -*- coding: utf-8 -*-
import math as m

p = [n*(3*n-1)//2 for n in range(0,1000)]
k = 0
sum_diff_pair = {}

while True:

	k += 1
	
	if k >= len(p):
		p += [n*(3*n-1)//2 for n in range(len(p),len(p)*2)]
	
	if p[k] in sum_diff_pair:
		print(sum_diff_pair[p[k]])
		break
	
	for j in range(1, k):
		p_l = p[k] + p[j]
		D = p[k] - p[j]
		l = int((m.sqrt(24*p_l+1)+1)/6)
		
		if p_l == l*(3*l-1)//2 and D in p:
			sum_diff_pair[p_l] = D

  어떤 숫자가 나왔고 그것이 답이었다. 하지만 그냥 넘어가기엔 뭔가 걸리는 부분이 있었다. 과연 최소 합이 최소 차를 보장하는가? 그렇지 않다. 50과 20의 합은 70이고 차는 30이다. 그리고 50과 50의 합은 100이고 차는 0이다. 과연 이 문제가 이와 같은 성질을 가지지 않는다는 보장이 있는가? 확신할 수 없었다.

  다음으로 생각한 것은 어떤 \(P_i\)가 \(D\)가 되기 위한 조건을 충족시키는지 확인하는 방법이었다. \(i\)를 순차적으로 늘려가면서 \(P_i\)를 만들 수 있는 \(k\)와 \(j\)가 존재하는지, 만약 그렇다면 \(P_l=P_k+P_j\)를 만족하는 \(l\)이 정수인지 확인하였다. 문제는 \(k\)의 범위였다. 최소값은 \(i+1\)이 명백하다. 하지만 최대값은 어디인가? 차의 최소값은 \(P_n-P_{n-1}=3n-2\)이므로 \(k\)의 최대값은 \((P_i+2)/3\)까지이다. 구현은 다음과 같다.

# -*- coding: utf-8 -*-
import math as m

i = 0

def p(n):
	return n*(3*n-1)//2
	
while True:

	i += 1
	p_i = p(i)
	print(i)
	for k in range(i+1, (p_i+2)//3+1):
		
		p_k = p(k)
		
		d = int((6*k - m.sqrt((1-6*k)**2 - 24*p_i) - 1)/6)
		p_j = p(k-d)
		
		tmp_p_i = p_k - p_j
		
		if tmp_p_i == p_i:
			
			p_l = p_k + p_j
			l = int((m.sqrt(24*p_l+1)+1)/6)
		
			if p_l == l*(3*l-1)//2:
				print(p_i)
				quit()

  이 방법으로 찾은 \(P_i\)는 최소임이 보장되며 답은 5482660이다. 하지만 큰 문제가 있는데, 그것은 이 방법은 너무 느리다는 것이다. 찾아야 하는 \(k\)의 범위는 \(O(i^2)\)이기 때문이다. 범위를 줄일 수 있는 방법을 깊이 고민했지만 결국 찾지 못했다.

  사람들은 이 문제를 어떻게 해결했는지 궁금하여 게시판에 들어가 보았다. 여러가지 답변이 있었지만 내가 원하는 답 즉, \(P_i=D\) 임을 보장하며 속도도 빠른 방법을 찾는 것은 쉽지 않았다. 그러던 중 한 프랑스인으로 추정되는 사람이 남긴 답변이 내 눈길을 사로잡았다. 그는 그가 제시한 방법 이전의 것은 모두 틀렸다고 단언하면서 다음과 같은 방법을 제시하였다.

  \(k=j+d\)라고 하면 \(D\)를 간단한 연산을 통하여 다음과 같이 쓸 수 있다.

$$D=P_k-P_j=P_{j+d}-P_j=3jd+P_d$$

  따라서 \(j=(D-P_d)/(3d)\)임을 알 수 있다. \(j\)는 정수이고 \(d\)의 범위는 \([1,i)\)이다. 이것으로 탐색해야 하는 범위가 획기적으로 줄어든다. 합 또한 \(P_l=j(3j-1)+D\)으로 간단하게 쓸 수 있다. 구현은 다음과 같다.

# -*- coding: utf-8 -*-
import math as m

i = 0

def p(n):
	return n*(3*n-1)//2

while True:
	i += 1
	D = p(i)
	
	for d in range(1, i):
		
		p_d = p(d)
		j, r = divmod(D - p_d, 3*d)
		
		if r == 0 and (m.sqrt(24*(j*(3*j-1)+D)+1)+1)%6 == 0:
			print(D)
			quit()

  답은 같지만 속도는 훨씬 더 빠르고 그 답이 옳다는 것 또한 보장된다. 이 방법을 사용하면 최소의 차 뿐만 아니라 이후의 차 또한 빠르게 구할 수 있다. 정말 간단한 방법이지만 나는 이것을 생각해 내지 못했고 그 사람이 작성한 방법을 이해하는 데에도 오랜 시간이 걸렸다. 그 사람은 이것을 어떻게 생각했을까? 내가 이것을 생각하기 위해서 필요했던 것은 무엇이었을까? 역시 세상은 너무 넓다는 생각이 든다.