ナルシシスト数とは
ウィキペディアによれば以下のような数をナルシシスト数と呼ぶそうです。
n 桁の自然数であって、その各桁の数の n 乗の和が、元の自然数に等しくなるような数をいう。例えば、$1^3 + 5^3 + 3^3 = 153$ であるから、153 はナルシシスト数である。
Narcissistic Number(Wolfram Mathworld)により詳しい説明と実際の数の列挙があります。ナルシシスト数は有限で理論上60桁以下しかないそうですが、実際の最大値は39桁までの88個だそうです。
ナルシシスト数を求める
ナルシシスト数を求めるコードはA005188のPythonのコードを参考にして書きました。
このコードのポイントは各桁の数の n 乗の和は順序に関係がないので「重複あり組み合わせ」のcombinations_with_replacementを使って計算した結果の数を、桁ごとにばらしてソートしたものが元の組み合わせと同じになるかという判断で行っています。
結果の表示ではナルシシスト数と元の組み合わせのタプルを表示しています。
from itertools import combinations_with_replacement
N = 10
for k in range(2,N+1):
pdl = [d**k for d in range(9+1)] # store d**k of 0..9
for cd in combinations_with_replacement(range(10), k):
pds = sum([pdl[pd] for pd in cd])
if tuple(int(d) for d in sorted(str(pds))) == cd:
print(k, pds, cd)
#3 370 (0, 3, 7)
#3 407 (0, 4, 7)
#3 153 (1, 3, 5)
#3 371 (1, 3, 7)
#4 8208 (0, 2, 8, 8)
#4 1634 (1, 3, 4, 6)
#4 9474 (4, 4, 7, 9)
#5 93084 (0, 3, 4, 8, 9)
#5 92727 (2, 2, 7, 7, 9)
#5 54748 (4, 4, 5, 7, 8)
#6 548834 (3, 4, 4, 5, 8, 8)
#7 9800817 (0, 0, 1, 7, 8, 8, 9)
#7 4210818 (0, 1, 1, 2, 4, 8, 8)
#7 1741725 (1, 1, 2, 4, 5, 7, 7)
#7 9926315 (1, 2, 3, 5, 6, 9, 9)
#8 24678050 (0, 0, 2, 4, 5, 6, 7, 8)
#8 24678051 (0, 1, 2, 4, 5, 6, 7, 8)
#8 88593477 (3, 4, 5, 7, 7, 8, 8, 9)
#9 146511208 (0, 1, 1, 1, 2, 4, 5, 6, 8)
#9 912985153 (1, 1, 2, 3, 5, 5, 8, 9, 9)
#9 472335975 (2, 3, 3, 4, 5, 5, 7, 7, 9)
#9 534494836 (3, 3, 4, 4, 4, 5, 6, 8, 9)
#10 4679307774 (0, 3, 4, 4, 6, 7, 7, 7, 7, 9)
n進数のナルシシスト数
ナルシシスト数はn進数でも考えることが出来ます。例えば3進数で"12"は$1^2+2^2=5=12_3$となりナルシシスト数です。Narcissistic Number(Wolfram Mathworld)でも9進数までがリストされています。
Pythonでもnumpyのbase_reprとint関数のbase指定をすると簡単にn進数（36進数まで）を扱うことが出来るので。上のプログラムを少し変えただけでn進数のナルシシスト数を求めることが出来ました。
from itertools import combinations_with_replacement
import numpy as np
for D in [3,12,16,19,24,30,36]: #Base
nar = []
for k in range(2,5+1):
pdl = [d**k for d in range(D+1)] # store d**k of d=0..D
for cd in combinations_with_replacement(range(D+1), k):
pds = sum([pdl[pd] for pd in cd])
if tuple(int(d,D) for d in sorted(np.base_repr(pds, D))) == cd:
#print(k, pds, np.base_repr(pds, D), cd)
nar.append(np.base_repr(pds, D))
print(f"D={D}: {nar[:11]}")
#D=3: ['12', '22', '122']
#D=12: ['25', 'A5', 'A83', '577', '668', '14765', '938A4']
#D=16: ['208', '5B0', '60B', '8C0', 'EA0', '173', '156', '5B1', '8C1', 'EA1', '248']
#D=19: ['9A', 'AA', '180', '6D0', '70D', '181', '6D1', '278', '58B', 'G8B', 'FDD3']
#D=24: ['J5G']
#D=30: ['28', 'S8', 'JBL', 'GEK', 'B4N9']
#D=36: ['8M0', 'GQ0', '8M1', 'GQ1', '48H', 'N7Q', 'L9Q', 'T40LU', '2L1C6', '6D1MM', '1HDBE']
このアルゴリズムはProject Euler, Problem 749: Near Power Sumsを解くのに役に立ちます。
（開発環境：Google Colab）