PyPyを試してみた @ Project Euler Problem 92

PyPyというのは高速なPythonのことです（テキトー）。
公式サイトに実行形式バイナリがあるのですぐに使える。
PyPyをvirtualenvの中で使うとか、pipやipythonを使うという話も、誰かが書いてくれている（他力本願）。
以下、普通のPythonをCPython（Cで実装されたPython，Classic-Python）と書きます。

今回は、CPythonでやるとどうも遅かったProject EulerのProblem 92に対して、PyPyを試してみた。

環境
OS：Windows 7 Enterprise 64-bit SP1
CPU：Intel Mobile Core 2 Duo P8400 @ 2.26GHz
メモリ：4.00 GB デュアル-Channel DDR3 @ 532 MHz

Project Euler
プログラミングで解くちょっと数学的な問題がたくさんあるサイト。他のプログラミングコンテストみたいな制限時間はないけど、公式のF&Qによると、実行して1分以内に答えが出なければ諦めてアルゴリズムを見なおしたほうが良いらしい。1分と言ってもCとPythonではかなり違ってくるが、おそらく問題としてはPythonでも上手く書けば1分を切れるように作ってあるのだろう。

Problem 92
どんな自然数でも各桁を二乗して足すことを繰り返せば 1 か 89 のどちらかになる（不思議だね）。10,000,000 未満の自然数のうち 89 になるものはいくつあるか？という問題。

要するに、すべての数をどちらのグループに属するのか判定しようって話。1 か 89 になるまで繰り返すと言うけど、各桁の自乗和が最大になるのは 9,999,999 のときの 9^2 * 7 = 567 なので、567 までの振り分けが出来れば、あとはそれぞれ一度だけ自乗和を求めることでどちらに属するか分かる。ということで次のようなコードを書いた。

count = 0
memo = [-1 for i in xrange(568)]
memo[1] = 0
memo[89] = 1
for i in xrange(2, 10000000):
    n = i
    while 1:
        n = sum([(ord(c)-48)*(ord(c)-48) for c in str(n)])
        if memo[n] != -1:
            break

    if memo[n]:
        count += 1
        if i < 568:
            memo[i] = 1
    elif i < 568:
        memo[i] = 0

print count

memo に 567 までの結果を入れて、あとはそれを使って判定している。8行目、

n = sum([(ord(c)-48)*(ord(c)-48) for c in str(n)])

これは

n = sum([int(c)**2 for c in str(n)])

と同じ意味で、若干の高速化を狙っている。というか ord() の方がけっこう速い。

結果
以下では何度か実行してかかっただいたいの時間を大まかに述べます。

上のコードをCPythonで実行すると、49秒ほどかかった。遅いねえ、アルゴリズムを改良したほうが良いんじゃない？だが断る！PyPyを使う。
そしてPyPyの結果は、6.9秒前後で7倍以上速い！

まとめ
PyPyはCPythonと比べるとだいぶ速い。場合によってはCよりも速いこともあるらしいけど、それはまた別の話。ちなみに上のコードをCにベタ移植すると、3.7秒くらいで解けた。
PyPyはJIT コンパイラなるもので実行時にコンパイルしているらしいけど、使用感は普通のPythonと何も変わらないところが凄い。面白かった記事

・流行りのJIT コンパイラは嫌いですか? — PyPy Advent Calendar 2011 v1.0 documentation
http://www.longsleeper.com/pypyadvent11

おまけ
PyPyはインライン展開もやってるらしいが、手動で高速化するなら、またまた8行目、

n = sum([(ord(c)-48)*(ord(c)-48) for c in str(n)])

この、いちいちリストを作って総和を求めるという操作を展開して、

tmp = 0
for c in str(n):
    tmp += (ord(c)-48)*(ord(c)-48)
n = tmp

と書き換えると、PyPyの結果は、4.3秒そこらまで速くなった。ただただ冗長で仕方がない。この変換をPyPyでやってくれたら良いと思う。

あと散々言われてることだけど、パイパイという名前は日本ではおあつらえ向きじゃないよなあ…

matsulibの日記

Ingredients as Code

PyPyを試してみた @ Project Euler Problem 92