콘텐츠로 건너뛰기
Home » 프로젝트 오일러 69

프로젝트 오일러 69

오일러의 피(phi)함수 φ(n)은 n보다 작거나 같으면서 n과 서로 소인 숫자의 개수를 나타냅니다. 예를 들어, 9보다 작거나 같으면서 9와 서로 소인 것은 1, 2, 4, 5, 7, 8이므로 φ(9)=6이 됩니다.

n서로 소인 수φ(n)n/φ(n)
2112
31, 221.5
41, 322
51, 2, 3, 441.25
61, 523
71, 2, 3, 4, 5, 661.1666…
81, 3, 5, 742
91, 2, 4, 5, 7, 861.5
101, 3, 7, 942.5

위에서 보듯이 n ≤ 10인 경우 n/φ(n)의 값은 n=6일때 최대가 됩니다.

n이 1,000,000 이하의 자연수일 때, n/φ(n)는 언제 최대가 됩니까?

https://euler.synap.co.kr/problem=69

접근

n/φ(n)에 대해서 n이 크면 클수록, φ(n)는 작으면 작을수록 n/φ(n)은 최대값이 다가가게 된다. 문제에서 주어진 표를 보면 φ(n)은 일정하게 늘어나는 것도 아니고 그렇다고 특정한 주기를 갖는 것도 아니다.

오일러의 피함수는 자신보다 작으면서 서로 소인 자연수의 개수를 나타낸다. 즉 피함수 값이 작다는 말은 서로 소인 수들의 개수가 작다는 것을 의미한다. 근 N 이하의 모든 자연수를 소인수분해 했을 때, N과 특정한 소인수를 공유하는 수가 많으면 많을수록 φ(N)은 작아질 확률이 크다.

오일러 피함수는 자신보다 작은 서로 소인 수의 개수를 의미한다. n/φ(n) 의 값이 최대가 되기 위해서는 φ(n)의 값이 n대비 작을수록 유리하다. 피함수의 값이 작다는 말은 서로 소인 수의 개수가 작다는 것을 의미한다. 예를 들어 169는 소인수 13만 가지고 있으므로, 13의 배수인 13, 26, 39, 52, 65, 78, 91, 104, 117, 130, 143, 156을 제외하고 156개의 자연수가 서로소이다. 따라서 φ(169) = 156 이다. 반면 210의 경우, φ(210) = 48 밖에 되지 않는데, 210은 2, 3, 5, 7 의 네 개의 소인수를 가지고 있기 때문이다.

즉 작고 다양한 소인수의 곱으로 구성된 수 일수록 그보다 작은 수와 서로 소가 될 가능성이 줄어든다. 따라서,  n/φ(n) 이 최대가 되는 백만 이하의 값은 연속된 소수의 곱이면서 백만을 넘지않는 수가 되어야 할 것이다. 이 때 연속된 수는 2부터 시작한다. 3부터 시작해서는 많은 짝수들과 서로 소가 되는 경우가 생기기 때문에 문제의 조건을 만족하는 최대값이 되기가 어려울 것이다.

연속된 소수의 곱이 1백만을 넘지 않는 최대값은 510510으로, 이는 2 * 3 * 5 * 7 * 11 * 13 * 17 이다. 그리고 이 때의 phi(n) 92160으로 n/φ(n) 의 값은 5.54 가량 된다.

from functools import reduce

product = lambda xs: reduce(lambda x, y: x * y, xs, 1)

def sieve(bound):
  s = [True] * (bound // 2)
  for i in range(3, int(bound**0.5) + 1, 2):
    if s[i // 2]:
      s[i * i // 2 :: i] = [False] * ((bound - i * i - 1) // (2 * i) + 1) 
  return [2] + [2 * i + 1 for i in range(1, bound // 2) if s[i]]

s = sieve(1_000_000)
i = 0
while True:
  if product(s[:i+1]) > 100_0000:
    print(product(s[:i]))
    break
  else:
    i += 1

참고로 phi 함수는 다음과 같이 구현할 수 있을 것인데, 최대 공약수 구하는 함수 자체가 느리기 때문에 순수하게 phi 함수를 구현하려했다가는 수십분이 걸릴 수 있을 것이다. (답이 되는 510510의 n / phi(n) 만 계산하는데도 0.5초 가량이 소요됨)

def gcd(a, b):
  if a < b:
    a, b = b, a
  while True:
    q, r = divmod(a, b)
    if r == 0:
      return b
    a, b = b, r

phi = lambda n: len([x for x in range(1, n) if gcd(n, x) == 1])
foo = lambda n: n / phi(n)