すこしふしぎ．

VR/HI系院生による技術ブログ．まったりいきましょ．（友人ズとブログリレー中．さぼったら焼肉おごらなきゃいけない）

【python】開き直り数を求める【便乗してみた】

python

こんばんは，1000chです． 先日ももえさんの記事で開き直り数というものを初めて知りました．

曰く，

３４３５という数は、１桁ずつに分けて、それぞれの数を自分と同じ数でべき乗して足し合わせると，元の数である３４３５になる。 　3435 ＝ 33 + 44 + 33 + 55 このような数を“開き直り数”と呼ぶことにする。ちょっと考えれば1(=11)もこの種の数であることがわかる。では１以上の整数で、この１と３４３５以外の開き直り数をすべて見つけていただきたい。ただし、ここでは00は0とする。

ほう．面白そう． ということでpython使って解いてみます．

開き直り数

開き直り数を一般化すれば，

ある数が与えられ，その各桁のベキ乗数の和がもとの数と等しくなる数

となります．ルール自体は単純ですね．

各桁で最も大きいのは9の9乗なので，開き直り数には上限があると思われます．そこでまずは

\( 9^{9} \geq 10^{n} - 1 \)

を解いて，開き直り数の最大を求めましょう．こんな感じですかね．

a = 9**9
b = 10
d = 1
print(a,b,d)
while a > b:
    a += 9**9
    b *= 10
    d += 1
    print(a,b,d)

実行結果

387420489 10 1
774840978 100 2
1162261467 1000 3
1549681956 10000 4
1937102445 100000 5
2324522934 1000000 6
2711943423 10000000 7
3099363912 100000000 8
3486784401 1000000000 9
3874204890 10000000000 10

というわけで，10の10乗以下の数で考えれば良さそうです．

目指せ１行

探索する最大値も見えたので，それでは実際にコーディングしてみましょう． ルール的に全探索のコードは簡単に書けそうです．

各桁のそれぞれについてベキ乗して和をとって，くらいであればリスト内包で行けるのでは？ ということで書いてみました．

ans = [x for x in range(10**10) if x == sum(map(lambda v:int(v)**int(v),list(str(x))))]

数値の各桁をlist化するところがちょっと回りくどいかも．いい書き方無いのかなぁ．

さて書けたのでとりあえず動かしてみます．

....が．

....終わりません．

一行にループが終わりません．試しに

import time
start = time.time()

ans = [x for x in range(10**5) if x == sum(map(lambda v:int(v)**int(v),list(str(x))))]

print(ans, time.time() - start)

くらいで動かすと

[1, 3435] 0.7693190574645996

と出力されました．3435が出力されているので，論理的には間違ってなさそうです． 実行時間的には0.7秒というところです．

ちょっと冷静に計算量を考えてみます．

\( 10^{10} \div 10^{5} = 10^{5} \)

\( 10^{5} \times 0.769 = 76900[sec] \)

\( 76900[sec] = 1281.6 [min] = 21.36[h] \)

なん...だと...

そりゃなかなか終わらない訳です．

ちょっとアルゴリズムを変える必要がありそうです．

重複組み合わせを使う

冷静に考えると，先ほどのアルゴリズムでは同じ数字の組み合わせを何度も検証しています．

例えば，123 / 132 / 213 / 231 / 312 / 321 などはすべて同じ和(1+4+27=32)になりますね．

先ほどは「ある数字」が「各桁べき乗数の和と等しいか」を検証していますが，

「桁の組み合わせ」が「開き直り数を満たすか」を検証した方が良さそうです．

例えば，

  • (1,2,3) > 1 + 4 + 27 = 32 > (3,2) ... ダメ
  • (3,3,4,5) > 27+27+256+3125=3435 > (3,4,3,5) ... ok

という感じです．いわゆる重複組み合わせというやつですね．

すこし計算量を考えてみましょう．0~9からnコ選ぶ重複組み合わせは

\( \frac{ ((10-1)+n)! } { (10-1)!n! } \)

で計算できるはず．試しにn=11で計算すると167960コ(10の5乗オーダー)の組み合わせがあるようです．

n=1~11で最大限考えてもせいぜい10の6乗オーダなので，10の10乗に比べれば相当時間短縮ができそうです．

ということで実際にコードを書いてみます．

import itertools 
hash = list(map(lambda x:x**x,[0,1,2,3,4,5,6,7,8,9]))
hash[0] = 0

for i in range(1,11):
    for j in list(itertools.combinations_with_replacement(range(10),i)):
        t = tuple(map(int,sorted(list(str(sum(list(map(lambda x:hash[x],j))))))))
        if t == j:
            print(sum(list(map(lambda x:hash[x],j))))

１行で書こうと思ったらさすがに長過ぎたので分けました． それでも7行目がカオスです． なんかlisp系見てるみたい．

重複組み合わせに関してはpythonのitertool内にそのものがあったので使ってみました． 仕様としてはイテレータを進めるごとに組み合わせのタプルがソートされた状態で返ってくるようです． であれば，そのタプルの要素をそれぞれベキ乗して和をとり，それを再び各桁のタプルを生成して比較すればよさそう．

てことで7行目では，

  • 桁数字のタプルからベキ乗した数値のリストを作る：list(map(lambda x:hash[x],j))
  • リスト中の数字の総和をとる：sum()
  • その和の各桁をリスト化する：list(str())
  • リストをソートする：sorted()
  • タプルに変換する：tuple(map(int,<sorted_list>))

という処理を行ってます．みづらいね！

こんな感じで生成されたタプルともとのタプルを比較，一致したら開き直り数ということになります．

実行結果

0
1
3435
438579088

実行時間は1.9941349029541016[sec]でした．すばらしい！

まとめ

  • ももえさんの記事に便乗した
  • 開き直り数

  • 計算量ってちゃんと考えないといけない

  • はてぶでtex表記つかえるようになった！(｀・ω・´) <- 今回一番言いたかったこと