Python
sympy

Python (SymPy) で完全数、友愛数、婚約数、双子素数を求める

はじめに

この記事は 博士の愛した数式 を愛した僕の備忘録です。作中に登場した数に関する概念を SymPy を使って計算してみます。

コード

準備

import sympy


def sum_divisors_except_itself(n):
  '''
  自分自身を除く正の約数の和を返す。
  '''
  return sum(sympy.divisors(n)[:-1])


def sum_divisors_except_one_and_itself(n):
  '''
  1 と自分自身を除く正の約数の和を返す。
  '''
  return sum(sympy.divisors(n)[1:-1])

完全数

完全数 (perfect number)とは、その数自身を除く約数の和がその数自身と等しい自然数のこと。

results = []
for n in range(1, 10000 + 1):
  if n == sum_divisors_except_itself(n):
    results.append(n)

results
# [6, 28, 496, 8128]

友愛数

友愛数 (amicable numbers) とは、異なる 2 つの自然数の組で、自分自身を除いた約数の和が互いに他方と等しくなるような数のこと。

results = []
for m in range(1, 10000 + 1):
  n = sum_divisors_except_itself(m)
  if n >= m:
    continue
  if sum_divisors_except_itself(n) == m:
    results.append((n, m))

results
# [(220, 284), (1184, 1210), (2620, 2924), (5020, 5564), (6232, 6368)]

婚約数

婚約数 (betrothed numbers) とは、異なる 2 つの自然数の組で、1 と自分自身を除いた約数の和が互いに他方と等しくなるような数のこと。

results = []
for m in range(1, 10000 + 1):
  n = sum_divisors_except_one_and_itself(m)
  if n >= m:
    continue
  if sum_divisors_except_one_and_itself(n) == m:
    results.append((n, m))

results
# [(48, 75), (140, 195), (1575, 1648), (1050, 1925), (2024, 2295), (5775, 6128)]

双子素数

双子素数 (twin prime) とは、差が 2 である 2 つの素数の組のこと。

results = []
for m in sympy.primerange(sympy.prime(2), 10000 + 1):
  n = sympy.prevprime(m)
  if m - n == 2:
    results.append((n, m))

results
# [(3, 5), (5, 7), (11, 13), (17, 19), (29, 31), (41, 43), (59, 61), (71, 73), (101, 103), ..., (9857, 9859), (9929, 9931)]