青木日記

■ たらいまわしべんち最終章

しつこくたらいをまわしてみる。

うーむ、原先生ファイナルバージョンがかなり速いですねえ。 第三引数だけ遅延評価するとは凄いアイデアです。

ただ、やっぱりそれでも引数をデカくすると Hugs に負けるんですよー。

~/c/tarai % time ruby yield-tarai.rb                       aamine@harmony
768
ruby yield-tarai.rb  2.20s user 0.01s system 101% cpu 2.187 total
　
~/c/tarai % time runhugs tarai.hs                          aamine@harmony
768
runhugs tarai.hs  2.13s user 0.01s system 101% cpu 2.118 total

これは (768,384,0) で試した場合で、引数を倍にするともっと差が広がります。 笹田さん版も同じですね。 やはりどうやっても完全な遅延評価をしている Haskell 版には敵わないようです。

■ Ruby / memoise

ちなみに笹田さんの多段ハッシュ版は次のように書いたほうが短くて速いです。

$memo = {}
　
def tarai( x, y, z )
  (($memo[x] ||= {})[y] ||= {})[z] ||= tarai0(x,y,z)
end
　
def tarai0(x, y, z)
  if x <= y
  then y
  else tarai(tarai(x-1,y,z),
             tarai(y-1,z,x),
             tarai(z-1,x,y))
  end
end
　
puts tarai(768, 384, 0)
　
　
~/c/tarai % time ruby mhash-tarai.rb     # 笹田さんのオリジナル
768
ruby mhash-tarai.rb  15.21s user 0.25s system 100% cpu 15.452 total
　
~/c/tarai % time ruby mhash-tarai2.rb    # 青木改良版
768
ruby mhash-tarai2.rb  10.41s user 0.08s system 100% cpu 10.487 total

また負の整数がないと考えてよいのなら 配列を使うとさらに高速化できます。

$memo = []
　
# x だけは -1 になることがあるので +1 しておく。
def tarai( x, y, z )
  (($memo[x+1] ||= [])[y] ||= [])[z] ||= tarai0(x,y,z)
end
　
def tarai0(x, y, z)
  if x <= y
  then y
  else tarai(tarai(x-1,y,z),
             tarai(y-1,z,x),
             tarai(z-1,x,y))
  end
end
　
puts tarai(768, 384, 0)
　
　
~/c/tarai % time ruby mhash-tarai3.rb                      aamine@harmony
768
ruby mhash-tarai3.rb  8.15s user 0.09s system 100% cpu 8.230 total

まだ「第三引数のみyield」には敵わないものの、かなり速くなりました。

■ C / memoise (1)

C の力任せに計算するやつが Ruby memo 版にすら負けているのが気にくわなかったので、 C でも多次元配列に計算結果をキャッシュするようにしてみました。

#include <stdio.h>
#include <stdlib.h>
　
#define K1 192
#define K2 96
#define K3 0
　
static int tarai(int, int, int);
static inline int tarai0(int, int, int);
static void initialize_cache_table(int);
　
#define max3(a,b,c) max(max(a,b),c)
#define max(x,y) ((x) > (y) ? (x) : (y))
　
int
main(int argc, char **argv)
{
    initialize_cache_table(max3(K1, K2, K3) + 2);
    printf("%d\n", tarai(K1, K2, K3));
    exit(0);
}
　
static int *cache;
static int table_width;
　
static void
initialize_cache_table(int width)
{
    int i, j, k;
　
    table_width = width;
    cache = malloc(width * width * width * sizeof(int));
    if (!cache) {
        fprintf(stderr, "malloc failed\n");
        exit(3);
    }
    for (i = 0; i < width; i++) {
        for (j = 0; j < width; j++) {
            for (k = 0; k < width; k++) {
                cache[(i * width * width) + (j * width) + k] = -2;
            }
        }
    }
}
　
static int
tarai(int x, int y, int z)
{
    int result;
    int key;
    int w = table_width;
　
    key = ((x+1)*w*w) + ((y+1)*w) + (z+1);
    result = cache[key];
    if (result != -2)
        return result;
    result = tarai0(x, y, z);
    cache[key] = result;
    return result;
}
　
static inline int
tarai0(int x, int y, int z)
{
    if (x <= y) return y;
    return tarai(tarai(x-1, y, z),
                 tarai(y-1, z, x),
                 tarai(z-1, x, y));
}

これもかなり速かったんですが、キャッシュの密度が粗すぎるため、 (768,384,0) では malloc() が失敗してしまいました。がくり。

■ C / memoise (2)

やはり多次元配列では芸がなさすぎました。 イリッフィベクタにしてみます。

#include <stdio.h>
#include <stdlib.h>
　
#define K1 768
#define K2 384
#define K3 0
　
static void *new_cache(void*);
static int tarai(int, int, int);
static inline int tarai0(int, int, int);
　
#define max3(a,b,c) max(max(a,b),c)
#define max(x,y) ((x) > (y) ? (x) : (y))
　
void **cache;
int cache_size;
　
int
main(int argc, char **argv)
{
    cache_size = max3(K1, K2, K3) + 2;
    cache = new_cache(0);
    printf("%d\n", tarai(K1, K2, K3));
    exit(0);
}
　
static void*
new_cache(void *val)
{
    void **p;
    int i;
　
    p = malloc(sizeof(void**) * cache_size);
    if (!p) {
        fprintf(stderr, "malloc failed\n");
        exit(3);
    }
    for (i = 0; i < cache_size; i++) {
        p[i] = val;
    }
    return p;
}
　
static int
tarai(int x, int y, int z)
{
    void **p1, **p2;
    int result;
　
    if (!(p1 = cache[x+1])) {
        cache[x+1] = p1 = new_cache(0);
    }
    if (!(p2 = p1[y])) {
        p1[y] = p2 = new_cache(-2);
    }
    if ((result = p2[z]) == -2) {
        p2[z] = result = tarai0(x,y,z);
    }
    return result;
}
　
static inline int
tarai0(int x, int y, int z)
{
    if (x <= y) return y;
    return tarai(tarai(x-1, y, z),
                 tarai(y-1, z, x),
                 tarai(z-1, x, y));
}

ちっ、こっちのほうが多次元配列を使うより短かかったか。

~/c/tarai % gcc -O6 -finline-functions -omemo-tarai memo-tarai2.c
memo-tarai2.c: In function `tarai':
memo-tarai2.c:54: warning: passing arg 1 of `new_cache' makes pointer from integer without a cast
memo-tarai2.c:56: warning: assignment makes integer from pointer without a cast
memo-tarai2.c:57: warning: assignment makes pointer from integer without a cast

黙れコンパイラ。俺が正義だ！

~/c/tarai % time ./memo-tarai                              aamine@harmony
768
./memo-tarai  0.57s user 1.68s system 101% cpu 2.227 total
　
~/c/tarai % time runhugs tarai.hs                          aamine@harmony
768
runhugs tarai.hs  2.11s user 0.02s system 100% cpu 2.122 total

おしい。もうちょっとなんですけどね。 まだキャッシュの密度が薄いみたいなので、 データ構造を工夫すればまだ速くできそうです。

ただ、こうやって速くしたバージョンが本当に正しい値を計算してるか わからないんですよね。Haskell 版にしても Ruby yield 版にしても、 プログラムが正しそうであることはちょっと見れば確かめられます。 でもここで書いた C 版のプログラムが本当に正しいかどうかは ちょっと見ただけではわかりません。 少なくとも引数に負の数を与えると SEGV するわけだし。

■ C / memoise (3)

ええい、ここまでやったら行くとこまで逝け！ st_table を利用した C / memoise version 3。

#include <stdio.h>
#include <stdlib.h>
#include "st.h"
　
#define K1 768
#define K2 384
#define K3 0
　
static int tarai(int, int, int);
static inline int tarai0(int, int, int);
　
st_table *cache;
　
int
main(int argc, char **argv)
{
    cache = st_init_numtable();
    printf("%d\n", tarai(K1, K2, K3));
    exit(0);
}
　
static int
tarai(int x, int y, int z)
{
    st_table *t1, *t2;
    int result;
　
    if (!st_lookup(cache, x, &t1)) {
        t1 = st_init_numtable();
        st_add_direct(cache, x, t1);
    }
    if (!st_lookup(t1, y, &t2)) {
        t2 = st_init_numtable();
        st_add_direct(t1, y, t2);
    }
    if (!st_lookup(t2, z, &result)) {
        result = tarai0(x,y,z);
        st_add_direct(t2, z, result);
    }
    return result;
}
　
static inline int
tarai0(int x, int y, int z)
{
    if (x <= y) return y;
    return tarai(tarai(x-1, y, z),
                 tarai(y-1, z, x),
                 tarai(z-1, x, y));
}
　
　
~/c/tarai % gcc -O6 -Wall -finline-functions memo-tarai3.c st.c alloc.c -omemo-tarai
memo-tarai3.c: In function `tarai':
memo-tarai3.c:28: warning: passing arg 3 of `st_lookup' from incompatible pointer type
memo-tarai3.c:30: warning: passing arg 3 of `st_add_direct' makes integer from pointer without a cast
memo-tarai3.c:32: warning: passing arg 3 of `st_lookup' from incompatible pointer type
memo-tarai3.c:34: warning: passing arg 3 of `st_add_direct' makes integer from pointer without a cast
memo-tarai3.c:36: warning: passing arg 3 of `st_lookup' from incompatible pointer type
　
~/c/tarai % time ./memo-tarai                              aamine@harmony
768
./memo-tarai  1.34s user 0.07s system 101% cpu 1.387 total

やたー

■ ちなみに

真の最速はコレみたいです。

~/c/tarai % ghc -O -o tarai tarai.hs                       aamine@harmony
~/c/tarai % time ./tarai                                   aamine@harmony
768
./tarai  0.01s user 0.00s system 579% cpu 0.002 total

■ to_int

見逃してましたが、4/13 の中田さんのやつが結構速いですね。 基本的には lambda を使うのと同じですけど、to_int を使って むりやりオブジェクト指向ぽくしてるのが特徴ですか。

# ちょっと修正したものを再掲
　
class Tarai
　
  include Comparable
　
  def initialize(x, y, z)
    @x, @y, @z = x, y, z
    @t = 0
  end
　
  def tarai_inspect(level = 3)
    return "?" if (level -= 1) < 0
    s = "Tarai(" <<
        [@x, @y, @z].map {|x| Tarai === x ? x.tarai_inspect(level) : x}.join(", ") <<
        ")"
    s << ("%#d" % @t) unless @t.zero?
    s
  end
　
  def +(val)
    @t += val
  end
　
  def -(val)
    @t -= val
  end
　
  def <=>(val)
    to_int() <=> val.to_int
  end
　
  def coerce(val)
    case val
    when Integer
      [val, to_int()]
    else
      super
    end
  end
　
  def to_int
    @t += if !@x
            0
          elsif @x <= @y
            @y
          else
            Tarai.new(Tarai.new(@x-1, @y, @z),
                      Tarai.new(@y-1, @z, @x),
                      Tarai.new(@z-1, @x, @y)).to_i
          end
    @x = nil
    @t
  end
　
  alias to_i to_int
　
  def to_s
    to_i.to_s
  end
　
end
　
　
t = Tarai.new(*ARGV.map {|n| n.to_i })
puts "#{t.inspect} = #{t}"

■ たらいまわしパック

全部まとめました。

• http://i.loveruby.net/archive/tarai.tar.gz

思う存分たらいまわしてください。

■ 文字列を整数に変換

めも

# Ruby
"15".to_i        # チェックが甘い。e.g. "a".to_i == 0
Integer("15")    # チェックが厳しい。 e.g. Integer("a") は例外
　
; Scheme
(string->number "15")
　
 -- Haskell
read "15" :: Int
　
 -- Haskell (2) HaXml から抜粋
atoi :: String -> Int
atoi = foldl (\x y-> x*10 + ord y - ord '0') 0

■ Scheme / delay force

遅ればせながら shiro さんとこの Wiki も見てきました。

おおっ、delay/force なんてので遅延評価できるんですね。 Scheme って凄いなー。

■ まだ続くらしい

うーむ、原先生最新バージョンは……これは…… デフォルト引数を使うことで第三引数の無駄な評価を さらに一段減らしたんですね。す、凄すぎ。 これで (768,384,0) までは Hugs に勝てるようになりました。