Rubyでバイナリデータに対するrindex検索の挙動でハマったので調べたことメモ

自分の手元の環境でこんなことが起きた。

$ ruby -v
ruby 3.1.2p20 (2022-04-12 revision 4491bb740a) [arm64-darwin21]
$ irb
irb(main):001:0> "\x01\x80\x00\x00".index("\x01")
=> 0
irb(main):002:0> "\x01\x80\x00\x00".rindex("\x01")
=> 1

\x010 番目にしかないのだから、 .index でも .rindex でも 0 が返ってくるはずではないの??

先に結論

バイナリデータを扱うときには必ずEncodingを ASCII-8BIT に指定しておくこと!

きっかけ

roo というgem (記事時点で 2.9.0)を使って、Excelファイルを開こうとした。

require 'roo'

p Roo::Excelx.new("hoge.xlsx")

これは問題ないが、 #initialize の引数は file_or_strem ということでファイルを open したものを渡しても良いはず、と

require 'roo'

File.open("hoge.xlsx") do |f|
  p Roo::Excelx.new(f)
end

ということをすると以下のような謎のエラーを吐いて落ちる。

/Users/sugyan/.rbenv/versions/3.1.2/lib/ruby/gems/3.1.0/gems/rubyzip-2.3.2/lib/zip/entry.rb:365:in `check_c_dir_entry_static_header_length': undefined method `bytesize' for nil:NilClass (NoMethodError)

      return if buf.bytesize == ::Zip::CDIR_ENTRY_STATIC_HEADER_LENGTH
                   ^^^^^^^^^

で困っていたので他の人に聞いてみたところ「自分のところではそんな問題は起きていない」と言われ。 試しにDockerでRuby環境作って実行してみると、確かに問題なく動く…何故? と調べはじめた。

String#rindex の謎挙動

自分の環境と問題なく動く環境とで比較しながら見てみると、roo が使っている rubyzip (記事時点で 2.3.2) の中でファイル内容を読み取った結果の値が違っていた。

module Zip
  class CentralDirectory
    include Enumerable

    END_OF_CDS             = 0x06054b50

    ...

    def get_e_o_c_d(buf) #:nodoc:
      sig_index = buf.rindex([END_OF_CDS].pack('V'))
      ...

buf の中から 0x06054b50 をpackしたもの つまり "PK\x05\x06" のデータを末尾から探すために .rindex() を呼んでいるのだが、手元の環境と 問題なく動く環境ではその値が 1 ズレている。。

何故?と思って色々試しているうちに冒頭の例のようなものに辿り着いた。

もう少し深く追う

少なくとも \x80 などを含まないASCIIだけのものであればおかしな挙動にはならない。

irb(main):001:0> "\x01\x7F\x00\x00".rindex("\x01")
=> 0

また、これが増えるとズレもどんどん大きくなる。

irb(main):001:0> "\x01\x80".rindex("\x01")
=> 0
irb(main):002:0> "\x01\x80\x80".rindex("\x01")
=> 1
irb(main):003:0> "\x01\x80\x80\x80".rindex("\x01")
=> 2
irb(main):004:0> "\x01\x80\x80\x80\x80".rindex("\x01")
=> 3

逆から走査する際に異常な文字が検出された際に読み飛ばし、その結果の走査数を文字列長から引いたものを返したことにより値がおかしくなる、というのが予想される。

ソースを読んでみる。 rstring 関連はこのへんのようだ。

#ifdef HAVE_MEMRCHR によって実装が分かれている。 memrchr というのがglibcに入っているもので、自分のmacOS環境では使えなくて Docker内のLinux環境などでは使える、これによって環境によって挙動が変わったりする、のかもしれない。

で、使えない環境の方での実装は。 対象が見つかるまで while loop の中で rb_enc_prev_char で前の文字を走査している、という感じのようだ。encodingの影響を受けそう…?

    while (s) {
        if (memcmp(s, t, slen) == 0) {
            return pos;
        }
        if (pos == 0) break;
        pos--;
        s = rb_enc_prev_char(sbeg, s, e, enc);
    }

Encodingと実行環境

Encodingが関係しているようだったので色々調べてみた。

Rubyではバイナリデータ(バイト列)も String で扱う。

バイナリの取扱い

Ruby の String は、文字の列を扱うためだけでなく、バイトの列を扱うためにも使われます。しかし、Ruby M17N には直接にバイナリを表すエンコーディングは存在しません。このため、バイナリを String で扱う際には、ASCII 互換オクテット列を意味する ASCII-8BIT を用います。これにより、ASCII 互換であるこの String は 7bit クリーンな文字列と比較・結合が可能となります。

https://docs.ruby-lang.org/ja/latest/doc/spec=2fm17n.html

ということで、バイナリの場合は本来は ASCII-8BIT のEncodingを持つ文字列として検索しなければならないのに、手元の環境では UTF-8 のままで検索しようとしていた、というところに「使い方の問題」があったようだ。

irb(main):001:0> "\x01\x80\x00\x00".encoding
=> #<Encoding:UTF-8>

Encodingが ASCII-8BIT になっていれば、rindex でも正しい値を得ることができそうだ。 全然知らなかったが String#bASCII-8BIT の複製を得ることもできるらしい。

https://docs.ruby-lang.org/ja/latest/method/String/i/b.html

irb(main):001:0> "\x01\x80\x00\x00".rindex("\x01")
=> 1
irb(main):002:0> "\x01\x80\x00\x00".force_encoding(Encoding::ASCII_8BIT).rindex("\x01")
=> 0
irb(main):003:0> "\x01\x80\x00\x00".b.rindex("\x01")
=> 0

なるほど〜。

ではこの Encoding は何で決まるのか?というのも上記リンクに書いてある。

リテラルエンコーディング

文字列リテラル正規表現リテラルそしてシンボルリテラルから生成されるオブジェクトのエンコーディングスクリプトエンコーディングになります。

またスクリプトエンコーディングが US-ASCII である場合、7bit クリーンではないバックスラッシュ記法で表記されたリテラルエンコーディングは ASCII-8BIT になります。

さらに Unicode エスケープ (\uXXXX) を含む場合、リテラルエンコーディングUTF-8 になります。

複雑。。。

とにかく通常はスクリプトエンコーディングがまず大事のようだ。

スクリプトエンコーディング

スクリプトエンコーディングとは Ruby スクリプトを書くのに使われているエンコーディングです。スクリプトエンコーディングは マジックコメントを用いて指定します。スクリプトエンコーディングには ASCII 互換エンコーディングを用いることができます。 ASCII 非互換のエンコーディングや、ダミーエンコーディングは用いることができません。

現在のスクリプトエンコーディング__ENCODING__ により取得することができます。

さらに、magic commentによってこのスクリプトエンコーディングを決定でき、それが無い場合の挙動も書かれている。

マジックコメントが指定されなかった場合、コマンド引数 -K, RUBYOPT およびファイルの shebang からスクリプトエンコーディングは以下のように決定されます。左が優先です。

magic comment(最優先) > -K > RUBYOPTの-K > shebang

上のどれもが指定されていない場合、通常のスクリプトなら UTF-8、-e や stdin から実行されたものなら locale がスクリプトエンコーディングになります。 -K オプションが複数指定されていた場合は、後のものが優先されます。

通常のスクリプト-e かどうか、でも変わったりするのね…。そして特に指定ない場合は最終的にはlocaleが使われる。 ので、 LC_CTYPE などでも挙動が変わり得るようだ。

$ ruby -e 's="\x01\x80\x00\x00"; p __ENCODING__, s.encoding, s.rindex("\x01")'
#<Encoding:UTF-8>
#<Encoding:UTF-8>
1
$ ruby -Kn -e 's="\x01\x80\x00\x00"; p __ENCODING__, s.encoding, s.rindex("\x01")'
#<Encoding:ASCII-8BIT>
#<Encoding:ASCII-8BIT>
0
$ RUBYOPT="-Kn" ruby -e 's="\x01\x80\x00\x00"; p __ENCODING__, s.encoding, s.rindex("\x01")'
#<Encoding:ASCII-8BIT>
#<Encoding:ASCII-8BIT>
0
$ LC_CTYPE=C ruby -e 's="\x01\x80\x00\x00"; p __ENCODING__, s.encoding, s.rindex("\x01")'
#<Encoding:US-ASCII>
#<Encoding:ASCII-8BIT>
0

という具合に、何もしないと UTF-8 文字列として扱ってしまうが、オプションや環境変数によって リテラルの Encoding を変えることもできる。

こうやって正しく ASCII-8BIT のEncodingを持つ文字列として検索すれば、 String#rindex の値がズレるということはなさそうだ。

つまり再現条件は

手元の環境で String#rindex の値がズレたのは2つの要因があって

  1. memrchr が使えない環境でビルドしたRubyで実行していて
  2. ASCII-8BIT でないEncodingの文字列に対して rindex をかけていた

ということになる。2つ揃っていないと起きないものと思われる。

Rooの問題

前述の、Rooでエラーが起きる件について考える。

そもそもEncodingが不明で渡ってくるものに対して安易に rindex を使うものではない、ということで rubyzip 側に落ち度がありそうではある。が 現在の master branch では リリースされている 2.3.2 とは大きく変わっていて、もう同じことは起きないのかもしれない。詳しくは追っていないので分からないが、今回のと関連するissueがあった。

既にCloseされているが、今回調べた Zip::CentralDirectory は直接使用するものではなく Zip::File を使ってください、ということのようだ。

つまりこの場合、 Roo 側の使い方が悪い、ということになりそう。 で、Rooの方も調べてみると ちゃんとそれに対応しようとしていると思われる修正があった。

これが取り込まれれば手元の環境でも問題なく動くようになるかな…?

もしくは現時点でも、使う側が渡す引数のEncodingを明示的に指定してやることで一応回避できそうではある。

require 'roo'

File.open("hoge.xlsx", "r:ASCII-8BIT") do |f|
  p Roo::Excelx.new(f)
end

Rubyのバグではないの?

しかし Encodingを正しく指定していなかったために起きているとしても、 String#indexString#rindex も検索した対象の開始indexを返すものであるはずなのだから、それがズレるのはおかしいのではないのか…?という気はする。

さらに memrchr が使えるか否かのビルド環境によってだけで結果が変わる(ちゃんと確かめてないけど おそらくそう…)、というのも気持ち悪い。

かなり昔から現在の実装になっているようだし 今からでは挙動を変えづらそうではある。 仕様といってしまえばそれまでだけど…。 このあたりはRuby開発陣の方々の見解も聞いてみたいところではあります。

3.2

ちなみに Ruby 3.2 からは String#byteindexString#byterindex が追加されるそうです。今回のようなバイナリデータに対する検索にピッタリ使えそうですね。

spherical linear interpolation(slerp)によるlatent spaceでのnoise補間

memo.sugyan.com

の記事を書いてから、先行事例の調査が足りていなかったなと反省。 Latent Seed の Gaussian noise 間での morphing はあんまりやっている人いないんじゃないかな、と書いたけど、検索してみると普通に居た。

Stable Diffusion が公開されるよりも前の話だった…。

そしてこの morphing を作るためのコードも gist で公開されている

読んでみると、2 つの noise の間の値を得る方法として slerp という関数が使われている。

        for i, t in enumerate(np.linspace(0, 1, num_steps)):
            init = slerp(float(t), init1, init2)

https://gist.github.com/karpathy/00103b0037c5aaea32fe1da1af553355#file-stablediffusionwalk-py-L179-L180

この定義は以下のようになっている。

def slerp(t, v0, v1, DOT_THRESHOLD=0.9995):
    """ helper function to spherically interpolate two arrays v1 v2 """

    if not isinstance(v0, np.ndarray):
        inputs_are_torch = True
        input_device = v0.device
        v0 = v0.cpu().numpy()
        v1 = v1.cpu().numpy()

    dot = np.sum(v0 * v1 / (np.linalg.norm(v0) * np.linalg.norm(v1)))
    if np.abs(dot) > DOT_THRESHOLD:
        v2 = (1 - t) * v0 + t * v1
    else:
        theta_0 = np.arccos(dot)
        sin_theta_0 = np.sin(theta_0)
        theta_t = theta_0 * t
        sin_theta_t = np.sin(theta_t)
        s0 = np.sin(theta_0 - theta_t) / sin_theta_0
        s1 = sin_theta_t / sin_theta_0
        v2 = s0 * v0 + s1 * v1

    if inputs_are_torch:
        v2 = torch.from_numpy(v2).to(input_device)

    return v2

https://gist.github.com/karpathy/00103b0037c5aaea32fe1da1af553355#file-stablediffusionwalk-py-L101-L125

不勉強で知らなかったが、これは "spherical linear interpolation"、日本語では「球面線形補間」と呼ばれるのかな、という手法で、3DCGの世界などではクォータニオン(Quaternion)の補間としてよく使われていたりする、ようだ。

en.wikipedia.org

書かれている通り、元々はクォータニオンの補間の文脈で紹介されていたが、次元数に関係なく適用することができるということらしい。 空間上での原点から2つの点へのベクトルを考え、その2つがなす角  \Omega をドット積とノルムを使って求めることができる。 始点から終点まで、 0 から  \Omega へと角度を線形に変化させながら2点を結ぶ円弧上を移動させていく、という感じだろうか。

 \displaystyle Slerp(p_0,p_1;t) = \frac{\sin{[(1-t)\Omega]}}{\sin{\Omega}}p_0 + \frac{\sin{[t\Omega]}}{\sin{\Omega}}p_1

画像生成モデルでの latent space における補間としては 2016年くらいには提案されていたようだ。 そういえば GAN で遊んでいたときにちょっとそういう話題を聞いたことがあったような気もする…。

arxiv.org

しかしこれによって Gaussian noise として分散を保ったまま変化させられる、ことになるのか…? 数学的なことはちょっとよく分からない。 まぁ球面上を角度だけ変えて動いていると考えればノルムが変化しないから大丈夫なのかな??くらいの感覚でしかない…。

実際にこれを使って補間していったときにどのような変化になるか、前回の記事で書いていた sqrt を使うものと上述の slerp を使うもので ランダムな2点を繋ぐ軌跡を見てみた。2次元、3次元 空間上でそれぞれプロットすると以下のような違いがあった。

間隔が異なるだけでなく、明らかに軌道も異なるところを通るようになるようだ。

実装

引用したコードでも良いが、自分でもPyTorch前提で Slerp を書いてみる。 NumPy への変換はしなくても計算できそうだった。

def myslerp(
    t: float, v0: torch.Tensor, v1: torch.Tensor, DOT_THRESHOLD: float = 0.9995
) -> torch.Tensor:
    u0 = v0 / v0.norm()
    u1 = v1 / v1.norm()
    dot = (u0 * u1).sum()
    if dot.abs() > DOT_THRESHOLD:
        return (1 - t) * v0 + t * v1
    omega = dot.acos()
    return (((1.0 - t) * omega).sin() * v0 + (t * omega).sin() * v1) / omega.sin()

前述の数式の定義通りに実装するだけではあるが、やはり特別な場合に注意する必要はある。2点の原点からの方向がまったく同じか、もしくは正反対の方向のとき、理論上では単位ベクトルのドット積は 1-1 になる。 が、float値の多次元配列で計算してみると 多少の誤差が生じることがある。

for _ in range(10):
    v0 = np.random.normal([1, 4, 64, 64])
    v1 = 1.0 * v0
    print((v0 * v1 / (np.linalg.norm(v0) * np.linalg.norm(v1))).sum())
0.9999999999999997
1.0
1.0
0.9999999999999998
0.9999999999999999
1.0
1.0
0.9999999999999998
1.0
1.0000000000000002
for _ in range(10):
    v0 = np.random.normal([1, 4, 64, 64])
    v1 = -1.0 * v0
    print((v0 * v1 / (np.linalg.norm(v0) * np.linalg.norm(v1))).sum())
-1.0000000000000004
-1.0
-0.9999999999999998
-1.0
-1.0
-1.0
-0.9999999999999999
-1.0
-1.0000000000000002
-1.0

slerp では このドット積の値から arccos() で角度を求める必要があるが、例えば np.arccos() は 引数が -1.01.0 の間に収まっていない場合に nan を返すことになってしまう。 また、後でこの角度の sin() で割る処理があるので、ドット積が 1.0 で返せていたとしてもその後の計算で今度は inf が出てきたりする。 ので、こういうケースでは球面ではなく普通に線形で繋いでしまった方が良いだろうね、ということで DOT_THRESHOLD というものが使われて処理を分岐させているようだ。


※追記

指摘いただいたが、まったく同じ向きのときはともかく 正反対の方向の場合は線形補間が正しいとは限らない。 例えば2次元において (1, 0)(-1, 0) を補間する場合は原点を通らず (0, 1) もしくは (0, -1) を通る円軌道になった方が良かったり 3次元の球でいうと北極と南極を繋ぐのは赤道上を通る球面軌道であって欲しい、など。 これは2点のどちらかを適度にブレさせれば算出できる場合もあるかもしれないし、上記のように明らかに通るべき点があればまず両点からそこにまず補間してやるなど、色々な考え方がありそう。これもまた用途や場合によると思われる。 今回の Stable Diffusion で使っている乱数においては正反対向きはまず有り得ないものと思ってよさそうだ。


この小数点誤差がどれくらいになるのかは、要素数だったり 計算の順番だったり numpy で計算するか torch で計算するか などによっても変わってくるようだ。 どれくらいの値になったら線形と見做すか、の閾値(上述では 0.9995)は場合によると思われる。 少なくとも Stable Diffusion で使う [4, 64, 64] の shape で torch.randn() で生成される乱数同士の場合は、通常はドット積は 0.1 にすら届くことがまずないので、極端な話ではそれくらい低くても問題なさそうではある。 morphing の2点間で同一の noise を使う可能性がない、という前提であれば分岐なくしても問題ないくらい。

比較

で、実際に slerp を使うように変更すると morphing はどう変化するのか。 以前に載せた anime girl での morphing の slerp 版を作ってみた。

並べて比較

ちょっと変化のタイミングが変わったかな…? という程度で 大きな差は感じられない。

考察

結局のところ、Stable Diffusion においては Gaussian noise である限り何らかの画像を生成できるようになっているし、指定の2点の間をどう遷移しようと違いはないのかもしれない。 ただ slerp は線形に角度変化していくという観点でも、その遷移の移動量としてはより安定した変化が期待できるかな、という気はする。 むしろ前回の記事の手法が適当に思い付きでやった割にはそこそこ近いものが出来ていたのすごかったのでは…? というくらい。

また、今回は Gaussian noise 間でということを考えていたが slerp の手法自体は別にどこでも使えるものだとは思うので知っておいて損はないはず。また StyleGAN とかで遊ぶことがあったらそこでも使ってみても良さそう。

あと、prompt の embedding 結果に対する morphing でも slerp を使うことはできるはず、だが どうなんだろう。実験してみていないけど、A から B を遷移するのに まったく無関係な C を通過することになってしまったり 余計な変化が増えてしまうだけのような気もするので、こちらは線形に最短距離で遷移して刻み幅だけ調整するくらいが良いのではないかな、と思っている。実際にはやはりどちらでも大して変わらないかもしれないし prompt 次第かもしれない。試行錯誤するときの選択肢として持っておくと良さそうではある。

Stable Diffusionでmorphing

ということで少し触って遊んでみたのでメモ。

Stable Diffusion をザックリ理解

先月公開された Stable Diffusion。

stability.ai

高精度で美しい画像を出力できる高性能なモデルながら、Google Colab などでも手軽に動かせるし、 Apple silicon でもそれなりに動かせる、というのが魅力だ。

中身については 以下の記事の "How does Stable Diffusion work?" 以降のところが分かりやすい。

huggingface.co

図をそのまま引用させていただくと

という仕組みになっていて、受け取る入力は "User Prompt" と "Latent Seed" のみ。 前者が「どのような画像を生成するか」を決めて、後者で「どんなバリエーションでその画像を生成するか」を決めるような感じ。

User Prompt は [77, 768] の空間にエンコードされて、これを使って [4, 64, 64] の Gaussian noise を scheduler によって繰り返し denoise していくことで目標の画像のための "latent image representations" を生成していく。最後はこれを VAE(Variational Auto-Encoder) で拡大していくことで最終的な [512, 512, 3] の画像を得る。 この途中の scheduler による sampling がアルゴリズムによって結果の質が変わってきたりするし、回数が少なすぎると denoise が足りなくて汚い出力になったりする、というわけだ。

ともかく、重要なのはこの 2 つの入力だけで出力が決まるということ、そして prompt の入力は [77, 768] に embedding されたものが使われる、ということ。 prompt の文字列を工夫していくのも良いが、そこから embedding されて渡すものを直接指定してしまっても良いわけだ。 また、 Latent Seed の noise の方も少しずつ変えていくことで少しだけ違う感じの出力を得たりすることができそうだ。

自分は以前に StyeleGAN で latent space を変化させて生成画像の morphing などをやっていたので、それと同じようなことをやってみることにした。

Prompt 間での interpolation, morphing

まずは 2 つの異なる prompt から生成される画像の間を補間して繋いでみる。

(当然ながら、京都と東京の町並みを補間するからといって愛知や静岡の風景が生成されたりはしない。)

scripts/txt2img.py の中にある、prompt から画像を生成する部分のメインはここにある。

github.com

    c = model.get_learned_conditioning(prompts)
    shape = [opt.C, opt.H // opt.f, opt.W // opt.f]
    samples_ddim, _ = sampler.sample(
        S=opt.ddim_steps,
        conditioning=c,
        batch_size=opt.n_samples,
        shape=shape,
        verbose=False,
        unconditional_guidance_scale=opt.scale,
        unconditional_conditioning=uc,
        eta=opt.ddim_eta,
        x_T=start_code
    )
    x_samples_ddim = model.decode_first_stage(samples_ddim)

model.get_learned_conditioning() で与えられた prompt 文字列から [N, 77, 768]Tensor に embedding されたものが得られる。 これを sampler (ここではデフォルトで使われる DDIMSampler を使用している) に与えて ddim_steps 回数の sampling を実行して [N, 4, 64, 64] の denoise された結果が得られる。これを model.decode_first_stage() に与えることで最終的な画像に使われる値が得られるようだ。

sampler.sample() には色々なパラメータがあるが、とにかく重要なのは conditioning (c) と x_T (start_code) だけ。これを変化させることで生成画像をコントロールしていく。

start_code の方はここでは固定した値を使うようにすることで、「何を描くか」だけを徐々に変化させていく様子を作れる。一度だけ乱数を生成してそれを繰り返し使うようにすると良い。

で、 c の方は 2 つの異なる prompt からそれぞれ embedding された値を取り出して、線形に変化させていく。

指定した cstart_code から生成画像だけを得るような関数を書いておくとやりやすい。 model のロード方法などについては割愛。

from contextlib import nullcontext

import numpy as np
import torch
from PIL import Image

from ldm.models.diffusion.ddim import DDIMSampler
from ldm.models.diffusion.ddpm import LatentDiffusion


def get_device() -> torch.device:
    ...


def load_model() -> LatentDiffusion:
    ...


model = load_model(...)


def generate(
    c: torch.Tensor, start_code: torch.Tensor, ddim_steps: int = 50
) -> Image:
    batch_size = 1
    device = get_device()
    precision_scope = torch.autocast if device.type == "cuda" else nullcontext
    with torch.no_grad():
        with precision_scope("cuda"):
            with model.ema_scope():
                uc = model.get_learned_conditioning(batch_size * [""])
                shape = [4, 64, 64]
                samples_ddim, _ = DDIMSampler(model).sample(
                    S=ddim_steps,
                    conditioning=c,
                    batch_size=batch_size,
                    shape=shape,
                    verbose=False,
                    unconditional_guidance_scale=7.5,
                    unconditional_conditioning=uc,
                    eta=0.0,
                    x_T=start_code,
                )

                x_samples_ddim = model.decode_first_stage(samples_ddim)
                x_samples_ddim = torch.clamp(
                    (x_samples_ddim + 1.0) / 2.0, min=0.0, max=1.0
                )
                image = (
                    255.0 * x_samples_ddim.cpu().permute(0, 2, 3, 1).numpy()[0]
                ).astype(np.uint8)
                return Image.fromarray(image)

ちなみに、 sampler.sample() された結果の方を線形に繋いで変化させていくという手法もあるのだけど、これはもはやどんな画像が生成されるかほぼ決定された後の値なので、 morphing しても単なる画像合成のような感じにしかならなくて面白くはない。

指定した cstart_code で画像を生成する準備ができたら、あとはその入力を作っていくだけ。

def morph_prompts(prompts: Tuple[str, str], steps: int) -> None:
    start_code = torch.randn([1, 4, 64, 64], device=get_device())
    c0 = model.get_learned_conditioning(prompts[0])
    c1 = model.get_learned_conditioning(prompts[1])
    for i in range(steps + 1):
        x = i / steps
        c = c0 * (1.0 - x) + c1 * x
        img = generate(c, start_code)
        img.save(f"morphing_{i:03d}.png")

これらを繋げてアニメーションさせれば、2 つの異なる prompt 間の morphing が出来上がる。

ただ、似ているものならまだあまり違和感ないが あまりに異なる 2 つを morphing させようとすると、急激に変化してしまって面白くない。

embedding された空間がどんなものかは未知だが、ともかく A と B の 2 点間には必ず「denoise された結果 A になるもの」と「denoise された結果 B になるもの」が分断される地点がどこかに存在してしまう。 それは A B の中心かもしれないし、少しズレたところかもしれないが、そのあたりで急激な変化が起こり得る。 ので、中点に近い位置は出来るだけ細かい step で刻んだ方がよりシームレスな morphing になりやすいように感じた。 ので、単純な線形に繋ぐのではなく双曲線関数で刻み幅を微妙に変えながら作ってみることにした。

    a = np.arccosh(5.0)
    for i in range(steps + 1):
        t = i / steps
        x = sinh(a * (t * 2.0 - 1.0)) / sinh(a) / 2.0 + 0.5
        c = c0 * (1.0 - x) + c1 * x
        ...

それでもやっぱり急激な変化は捉えられないことが多々あるけれども…。

Seed 間での interpolation, morphing

今度は、同一の prompt で異なる Latent Seed を使用した 2つの画像間での morphing。

prompt の方は固定して、 torch.randn() で生成していた Gaussian noise の方を徐々に変えていく。生成する画像の「お題」は一緒だが、違うバリエーションのものになっていく、という morphing。

prompt のときと同じように変化させていけば良いだけ、と思ったが そうはいかない。実際やってみると中間点あたりはボヤけた画像になってしまうようだ。 最初 何故だろう…?と思ったが どうやらこの noise は "Gaussian noise" であることが重要で、標準正規分布として  {\mu} = 0, {\sigma}^2 = 1 になっていなければならない、ということらしい。 単純に v0 * (1.0 - x) + v1 * x のように単純な線形結合で変化させていくと、中心に近づくにつれてその標準偏差は小さくなってしまう。

それを防ぐために、足し合わせる前にそれぞれの倍率の sqrt をとるようにすると、合成された noise は標準偏差を保持したまま遷移することができそうだ。

すると今度は 0 付近と 1 付近で急激な変化が起こりやすそうなので、 prompt morphing のときのように刻み幅を調節する。

def morph_noises(prompt: str, steps: int) -> None:
    c = model.get_learned_conditioning([prompt])
    n0 = torch.randn([1, 4, 64, 64], device=get_device())
    n1 = torch.randn([1, 4, 64, 64], device=get_device())
    for i in range(steps + 1):
        t = i / steps
        x = 2.0 * t**2 if t < 0.5 else 1.0 - 2.0 * (1.0 - t) ** 2
        start_code = n0 * math.sqrt(1.0 - x) + n1 * math.sqrt(x)
        img = generate(c, start_code)
        img.save(f"morphing_{i:03d}.png")

これで、同じお題(prompt)に対して複数のバリエーションで描かれたものを連続的に変化させていくことができる。 好みの画像を出力する seed を幾つかピックアップして繋いでみたりするとより好みのものが見つかるかもしれないし、意外とブレンドされたものは好みではないものになるかもしれない。


※追記

memo.sugyan.com


まとめ

以上の2つができれば、その応用として prompt と noise を同時に変化させていったり、交互に変化させていったり、数回変化させた後にまた元の画像に戻ってきたり、といったものも作っていける。

雑に試行錯誤しながら実行できるように Google Colab でscriptを書いていたけど、ちょっと整理して公開する予定(需要あるかどうか分からないけど)。 prompt 間の morphing はやっている人結構いるけど、noise 間のものはまだあんまり見かけないような気はする?

※追記: GitHub - sugyan/stable-diffusion-morphing でとりあえず公開しておきました。多分動くはず?

Slackの家庭内日記をはてなブログに同期していく

夫婦で日記の習慣を始めたのは3年ほど前。簡単な記録でもつけておくと後で色々見返せて良いよね、と。

とはいえなかなか習慣をつけるのは難しいなと思って

  • ブログサービスは最も日付と関連づけて文書が保存されるし良い
  • しかしブログ投稿、となるとちょっと気合いが要る
  • SNS感覚で気軽に投稿できると良い
  • が、下手なSNS誤爆したりすると大変

といったものがあり、色々試行錯誤もした結果 家庭内Slackで日記用channelを作り、1日1投稿ずつ書き込む、という運用にしたのだった。

  • 誤爆の心配は少ない (少なくとも妻は家庭以外でSlackを使わないし)
  • 専用channelで投げていくだけなので他の話題と混ざったりもしない
  • 気軽に書けて編集なども簡単
  • 検索もそれなりにできて過去の記録を見返しやすい

ということで気に入っていた。

ところが90日以上前のものが見られなくなるということで移行を検討した。 べつに課金すれば良いだけではあるのだけど、せっかくならバックアップも兼ねて別のブログに同期しておいたら良いのでは、と。

ので、過去の日記投稿たちをSlackから取得してブログ記事として投稿していく、というものを作った。

なんとなく書きやすいのはPythonかな?と思ってPythonで。ライブラリも充実していて数十行程度でそれぞれのロジックが書けてラクだったので正解だったと思う。

APIでの取得と投稿

SlackのAPIからchannelの履歴を遡り取得していく。 text要素だけ取ると絵文字が :bow: などで表示されてしまうので rich_text_section から復元していってSlackに書かれた通りの文字列でとれるようにした。

同期先のブログははてなブログにした(自分が使っていて多少は勝手が分かるので)。自分と妻のアカウントだけが見られる限定公開にして、そこに AtomPubAPIで投稿していく。 updated を指定すれば日時を指定できるので過去のものも過去の日記として記録されていく。

同期をするタイミングは確定していないけど、複数回実行するし 一度投稿したものも編集する可能性はゼロではないので、既に記事があるものは PUT で更新、まだ無いものは POST で新規投稿する。が 投稿しようとしているものが既に記事あるかどうかを調べるのが難しく、ドキュメントには entry_id がepochと書いてあるけどそうでもない… まぁ投稿日時が一致するようにするから updated に使う日時と同じ epoch の記事が既にあるかどうかを調べれば良いか、ということでそうした。

はてなブログの制限

いざ過去の日記をどんどん記事として投稿していこうとしたら50日ぶんくらいで止まってしまった。 知らなかったけどどうやらはてなブログには「1日100件まで」という上限があるらしい…。 まぁ3年分なら毎日コツコツやっていけば 365 * 3 * 2人 = 2190件 だから約21日で終わるはず。ギリギリ8月中には。。

もっと急ぎで大量の記事を移したかったら インポートできる形式のファイル を作るとかすることになるんだろうか。

定期実行

過去のはボチボチやっていくとして、今後これから書くものをどうするか。 すっかり3年の習慣になっているので今からブログに全面移行するつもりもなくて、日記は今まで通りにSlackで書き続けて、定期的に直近の日記をブログに同期していけば良い、ということにした。

Cloud Functionsにdeployして Cloud Schedulerから定期実行して Pub/Sub経由で自動で定期実行してくれるように設定。これであとは勝手に同期されるようになってくれるはず。

書くときはSlackに書き続け、90日以上前のものは調べたり検索したりできなくなるけど そのときはブログの方で検索すればもっと昔まで遡ることができる、という運用。

まとめ

まだ完全に同期に成功しているわけではないけど、この運用でしばらくやってみる予定。そもそも家庭内の日記がどこまで続くか…という問題もあるけど。まぁ夫婦ともに3年くらい続いているわけだし しばらくは続くんじゃないかな? 子育てで大変な日々も数年後には笑って読み返せるといいな

Rustで将棋棋譜変換ライブラリを作った

というわけで作った記録

将棋棋譜の形式色々

プロ対局、コンピュータ将棋、対戦アプリ、…と様々な将棋対局についてのデータが世の中には存在しているが、その形式は様々だったりする。

KIF

将棋ソフト「柿木将棋」 で使われている形式だそうで、現在でも多くの将棋ソフトなどで使われているようだ。プロ対局の中継サイトなどでもこの形式が使われている。

例:

# ---- Kifu for Windows95 V3.53 棋譜ファイル ----
開始日時:1999/07/15(木) 19:07:12
終了日時:1999/07/15(木) 19:07:17
手合割:平手
先手:先手の対局者名
後手:後手の対局者名
手数----指手---------消費時間-- # この行は、なくてもいい
1 7六歩(77) ( 0:16/00:00:16)
2 3四歩(33) ( 0:00/00:00:00)
3 中断 ( 0:03/ 0:00:19)

初期盤面や指し手など日本語で記述されていて、人間が読むぶんには直感的に分かりやすい。

が、コンピュータで解析するにはなかなか厄介…。 書式の定義については以下のリンクだけが存在している状態のようで、これを頼りに解析していくしかない。しかし「省略できる」「任意のものを追加できる」など なかなか自由な項目が多くて対応するのは大変。

優れている点としては、上記リンクでは言及されていないが「指し手の変化・分岐が表現できる」、という点で、枝分かれしたものも含めてすべて記録として残せるようになっているようだ。自分は柿木将棋などのシリーズ製品を持っていないけど、おそらくこれらがちゃんと分岐を含む棋譜も読み取って再現できるようになっていると思われる。

CSA

これは Computer Shogi Association で定めている形式のようで、KIFと比較すると圧倒的にコンピュータで解析しやすい形になっている。

例:

'----------棋譜ファイルの例"example.csa"-----------------
'バージョン
V2.2
'対局者名
N+NAKAHARA
N-YONENAGA
'棋譜情報
'棋戦名
$EVENT:13th World Computer Shogi Championship
'対局場所
$SITE:KAZUSA ARC
'開始日時
$START_TIME:2003/05/03 10:30:00
'終了日時
$END_TIME:2003/05/03 11:11:05
'持ち時間:25分、切れ負け
$TIME_LIMIT:00:25+00
'戦型:矢倉
$OPENING:YAGURA
'平手の局面
P1-KY-KE-GI-KI-OU-KI-GI-KE-KY
P2 * -HI *  *  *  *  * -KA * 
P3-FU-FU-FU-FU-FU-FU-FU-FU-FU
P4 *  *  *  *  *  *  *  *  * 
P5 *  *  *  *  *  *  *  *  * 
P6 *  *  *  *  *  *  *  *  * 
P7+FU+FU+FU+FU+FU+FU+FU+FU+FU
P8 * +KA *  *  *  *  * +HI * 
P9+KY+KE+GI+KI+OU+KI+GI+KE+KY
'先手番
+
'指し手と消費時間
+2726FU
T12
-3334FU
T6
%CHUDAN
'---------------------------------------------------------

主にASCII文字で表現されるようになっていて 駒種は2文字、場所を示すのも2桁の数字、といった一定の決まりがあるので解析しやすい。 これと同様の表現を用いたCSAサーバプロトコルがあり、世界コンピュータ将棋選手権などではこれが使われているらしい。

USI/SFEN

これは棋譜形式とはちょっと違うかもしれないけど、 将棋GUIソフト「将棋所」 で使われて広く普及しているUSIという通信プロトコルがある。 GUIと思考エンジンの間の通信をするために局面情報や指し手の情報を SFEN(Shogi Forsyth-Edwards Notation)という形式で表現していて、これによって棋譜を再現できる。

例:

position sfen lnsgkgsnl/9/ppppppppp/9/9/9/PPPPPPPPP/1B5R1/LNSGKGSNL w - 1 moves 5a6b 7g7f 3a3b

文字表記

そもそもコンピュータで扱われる前は文字情報だけで記録が残っているはずで、これは指し手をそれぞれ数字と駒種と相対位置や動作で表現している。KIF形式はこれらの表現をカスタマイズした形式といえる。

76歩
52銀右上成

その他

その他にも幾つかの形式があるようで、どれがどのくらい普及しているのかよく分からない…。 おそらく多くはKIFかCSAが使われていると思われる。CSAはむしろコンピュータ将棋の人にしか使われていないか…?

GUIソフトや将棋ソフトの開発者はそれぞれこれらの形式を読み書きできるように対応していくしかない、という状況のようだ…

JSON棋譜フォーマット(JKF, json-kifu-format)

上述した棋譜形式はそれぞれ読みやすさや表現できる情報が違っているわけだが、実際のところ指し手には多くの情報が必要で、単独の指し手の表記だけを取り上げたときに「どの駒が」「どこから動いて」「どうなったか」などの情報をすべて知ることは出来ない。直前の指し手や現在の盤面の情報がないと分からないものが多かったりする。

そこで、json-kifu-formatという形式が提案されて公開されている。

github.com

READMEに書かれている通り、各指し手について必要な情報をすべて含めたJSON形式で表現し、また分岐にも対応している。

例:

{
  "header": {
    "先手": "na2hiro",
    "後手": "うひょ"
  },
  "moves": [
    {},
    {"move":{"from":{"x":7,"y":7},"to":{"x":7,"y":6},"color":0,"piece":"FU"}},
    {"move":{"from":{"x":3,"y":3},"to":{"x":3,"y":4},"color":1,"piece":"FU"}},
    {"move":{"from":{"x":8,"y":8},"to":{"x":2,"y":2},"color":0,"piece":"KA","capture":"KA","promote":false}},
    {"move":{"from":{"x":3,"y":1},"to":{"x":2,"y":2},"color":1,"piece":"GI","capture":"KA","same":true}},
    {"move":{"to":{"x":4,"y":5},"color":0,"piece":"KA"}},

    {"special": "CHUDAN"}
  ]
}

さらにこの形式の提案だけでなく、KIFやCSAからJKF形式に変換するためのparserと、JKFを再生するためのプレイヤーも同梱している。

実装としては それぞれの棋譜形式を一度parseしたのち、初期局面から一手ずつ実際に進めて足りない情報を補っていくことで最終的なJKFを生成しているようだ。

自作parser, converter

上記 json-kifu-parser はJavaScriptで実装されているが、自分は最近はRustで将棋関連のものを作っていてRustでも同様のものが欲しいと思った。

具体的には、自作の詰将棋ソルバ に色々な問題を食わせてみようと思ったときにKIFファイルが多かったり 自分で書くにはCSAで書いたりで色々な入力に対応したいというのがあった。

CSAのparseには csa crateが既に存在していて問題ないが、KIFは自力でやるしかなく、いちおう自前である程度のparseできるところまで作っていたが、折角だから分岐もちゃんと対応してJKFに互換のデータを作れるParserを作ってみよう、という気持ちになった。 ので作ることにした。

github.com

KIFのparseには csa と同様に nom というparser combinators libraryを使って頑張って実装した。optmany0 が多くてなかなか大変だったが、どうにか様々なケースに対応できたと思う。 パフォーマンス的にちょっと問題があるかもしれないけど… 大量の棋譜を食わせて高速に処理をしなければならない、という場面でない限りは大丈夫か。

分岐は一度最後まで読み込んでstackに積んでおいてから後で forks に詰め替えていく必要があってなかなか頭を使った。

CSAのparseはcsa crateに任せてparse結果からJKFに変換する処理だけを書いて作ったが、csa crateは指し手のコメントなどを捨ててしまっていて復元できなかったりもあるので、こちらも自分で実装しなおしても良いのかもしれない。

それぞれparseして結果を JsonKifuFormat に詰めたら、足りない情報を補うために normalize() を呼ぶ。 この中で実際に局面を動かして確かめる必要があるので内部で shogi_core crateを使った。forks再帰的に処理していく。

また、 relative 情報は他の同種駒と比較して相対的なものを調べる必要があり、計算が大変めんどう…なので shogi_official_kifu crateで表記結果から取り出して詰めることでお茶を濁した。本来ならこの棋譜表記を作成するための情報なので本末転倒な感もあるが…。そもそも relative はそういった棋譜表記をするときくらいにしか使わない(棋譜変換などでは必要ない)ので計算もoptionalで良いのかもしれない。

一度 JsonKifuFormat が作れれば、そこからKIF形式やCSA形式に書き出すのは比較的容易にできる。(KIFの分岐はまだ対応できていないが…)。 また shogi_core::Position にも変換できるし shogi_core::ToUsi traitを利用してUSIで書き出すのもやりやすい。

ので、もし対応する棋譜形式を増やそうと思ったらそれぞれ JsonKifuFormat と相互変換するためのparser/converterだけを用意していけば良い、ということになる。 記事冒頭に載せた図の通りで、変換の中継のための形としてとても良いデータ型になっていると思う。

クレート公開

せっかく便利ライブラリとして実装できたので、(需要はまったく無いかもしれないが)crates.ioに公開してみることにした。 できるだけドキュメントをちゃんと書いて、初の cargo publish

docs.rs

Rustで将棋棋譜を扱いたい方が居たら使ってみていただけると嬉しいです。

WASM化してブラウザ上でclient sideだけで棋譜変換できるWebApp、なら多少は使いどころあるかな…?

SIMDによる将棋Bitboard計算の高速化

ということでSIMDでの高速化のメモ。

SIMDとは

ja.wikipedia.org

の通り、複数のデータを1命令で同時に演算する、というもの。 将棋Bitboardは81マスのデータを表現するのに64bitでは足りないので、一般的に [u64; 2] など複数の要素を使うことになる。 このBitboard同士の論理演算を普通に実装しようとすると、

struct Bitboard([u64; 2]);

impl std::ops::BitAnd for Bitboard {
    type Output = Self;

    fn bitand(self, rhs: Self) -> Self::Output {
        Self([self.0[0] & rhs.0[0], self.0[1] & rhs.0[1]])
    }

}

のように各要素同士をそれぞれ演算した結果を新たに格納する、という操作が必要になる。 こういった場面でSIMDを利用することで、u64同士の2つの演算を1命令で同時に行うことができる。

2命令が1命令に減る程度ではあるが、合法手生成のための駒の利きの計算などで多くの論理演算を行っているので、計算の高速化のためには重要なものになる。

実装

defaultの実装にしても良いのかもしれないが、とりあえずは features として定義した。

[features]
simd = []

cargo build --features simd など指定したときのみ有効になり、これが指定されなかったり 対応する実装が存在していない場合は、SIMDを使わない shogi_core::Bitboard を使用する実装にfallbackするようにしている。

x86_64

まずはx86_64のものを実装した。128bitレジスタを利用するSSE命令をメインに、一部でAVX2を使用して256bitレジスタでの計算も行っている。

use std::arch::x86_64;

pub(crate) struct Bitboard(x86_64::__m128i);

SIMD Bitboard by sugyan · Pull Request #17 · sugyan/yasai · GitHub

基本演算

前述したような論理演算については x86_64::_mm_and_si128, x86_64::_mm_or_si128, x86_64::_mm_xor_si128 でそれぞれ128bitの計算を1命令で出来るのでそれを使う。

同値判定、ゼロ値判定には x86_64::_mm_test_all_zerosi32 の値を返してくれるのでこれで判定できる。

飛び利き計算

前回の記事 でLeading/Trailing Zerosを使う方法について書いたが、このLeading/Trailing Zerosを求めるためには一度SIMDレジスタから何らかの形で値を取り出す必要があり、SIMDとはあまり相性が良くないようだった。ので極力使わずに論理演算のみで解決する方向を模索した。

まずは一番簡単な後手香車の利き。decrementして xor とる、という方法。 1〜7筋か/8,9筋か、で対応する u64 を取り出して -1 するのが一番命令数少なくて済みそうではあるが、分岐をなくして論理演算だけで計算するようにしてみた。

 impl Bitboard {

    ...

   fn sliding_positive_consecutive(&self, mask: &Self) -> Self {
        unsafe {
            let and = x86_64::_mm_and_si128(self.0, mask.0);
            let all = x86_64::_mm_cmpeq_epi64(self.0, self.0);
            let add = x86_64::_mm_add_epi64(and, all);
            let xor = x86_64::_mm_xor_si128(and, add);
            Self(x86_64::_mm_and_si128(xor, mask.0))
        }
    }
 }

_mm_cmpeq_epi64 に同値を渡しすべて1になっているものを作り、それを足すことで両方を同時に1減じる。 ここではmaskがある筋の連続的なbitしか渡ってこない前提であり、最終的にこのmaskとの論理積をとるので関係ない部分まで変化してしまっても問題ない。

先手香車はQugiyの手法が最もシンプルで良さそうだったのでそのまま使わせていただいた。

impl Bitboard {

    ...

    fn sliding_negative_consecutive(&self, mask: &Self) -> Self {
        unsafe {
            let m = x86_64::_mm_and_si128(self.0, mask.0);
            let m = x86_64::_mm_or_si128(m, x86_64::_mm_srli_epi64::<1>(m));
            let m = x86_64::_mm_or_si128(m, x86_64::_mm_srli_epi64::<2>(m));
            let m = x86_64::_mm_or_si128(m, x86_64::_mm_srli_epi64::<4>(m));
            let m = x86_64::_mm_srli_epi64::<1>(m);
            Self(x86_64::_mm_andnot_si128(m, mask.0))
        }
    }
}

飛車・角行の飛び利きは、4方向をそれぞれ(Squareのindex的に)正方向のものと負方向のもので分けて2つずつ計算する。 本当は飛車の縦方向は香車の利きを再利用するのが良いのだとは思うけど、両駒の利きを同じinterfaceで計算できるようにするためにそれはあえてしなかった。

正方向は前述のdecrementしてxor取るものを、ちゃんと桁借りを考慮して128bit全体のdecrementになるようにした上で、AVX2の256bitレジスタを使って2方向分を同時に計算する、という形に。

impl Bitboard {

    ...

    fn sliding_positives(&self, masks: &[Self; 2]) -> Self {
        unsafe {
            let self256 = x86_64::_mm256_broadcastsi128_si256(self.0);
            let mask256 = x86_64::_mm256_set_m128i(masks[0].0, masks[1].0);
            let masked = x86_64::_mm256_and_si256(self256, mask256);
            // decrement masked 256
            let all = x86_64::_mm256_cmpeq_epi64(self256, self256);
            let add = x86_64::_mm256_add_epi64(masked, all);
            let cmp = x86_64::_mm256_cmpeq_epi64(add, all);
            let shl = x86_64::_mm256_slli_si256::<8>(x86_64::_mm256_xor_si256(cmp, all));
            let dec = x86_64::_mm256_sub_epi64(add, shl);
            // (masked ^ masked.decrement()) & mask
            let xor = x86_64::_mm256_xor_si256(masked, dec);
            let ret = x86_64::_mm256_and_si256(xor, mask256);
            Self(x86_64::_mm_or_si128(
                x86_64::_mm256_castsi256_si128(ret),
                x86_64::_mm256_extracti128_si256::<1>(ret),
            ))
        }
    }
}

負方向はおそらくQugiyのように swap_bytes して unpack して桁借りを計算して… というのが良いのだと思いつつも、前述した leading_zeros を使うのが気に入っていたのでこれを応用することにした。

SIMDレジスタ上での leading_zeros の正確な値は求めづらいが、「leftmost set bitが含まれる8bit」がどこかは見つけやすい。 _mm256_movemask_epi8 を使って、各8bitにおける最上位の値を集めた u32 を得ることができる。 ので、 _mm256_cmpeq_epi8_mm256_setzero_si256 と比較した結果に対してこのmovemaskをかけてやれば、「各8bitにおいて値が 0 ではなかったもの」をbit列で得ることが出来る。あとは16bitずつ上位と下位で切り出して反転させてから leading_zeros をとれば、どの8bitに leftmost set bit があったかを知ることができる。 ということでその8bit以降を全部塗り潰すmaskだけ用意しておく。

あとは先手香車のように >>(shift right) と or を3回だけ繰り返すことで、各8bit内で leftmost set bit 以降をすべて 1 にすることは出来るので、それと上述のmaskの論理和をとれば「leftmost set bit以降がすべて 1 になっているもの」を作り出すことはできる。

const MASKED_VALUES: [(i64, i64); 16] = [
    (0, 0),
    (0x0000_0000_0000_00ff, 0),
    (0x0000_0000_0000_ffff, 0),
    (0x0000_0000_00ff_ffff, 0),
    (0x0000_0000_ffff_ffff, 0),
    (0x0000_00ff_ffff_ffff, 0),
    (0x0000_ffff_ffff_ffff, 0),
    (0x00ff_ffff_ffff_ffff, 0),
    (-1, 0),
    (-1, 0x0000_0000_0000_00ff),
    (-1, 0x0000_0000_0000_ffff),
    (-1, 0x0000_0000_00ff_ffff),
    (-1, 0x0000_0000_ffff_ffff),
    (-1, 0x0000_00ff_ffff_ffff),
    (-1, 0x0000_ffff_ffff_ffff),
    (-1, 0x00ff_ffff_ffff_ffff),
];

impl Bitboard {

    ...

    fn sliding_negatives(&self, masks: &[Self; 2]) -> Self {
        unsafe {
            let self256 = x86_64::_mm256_broadcastsi128_si256(self.0);
            let mask256 = x86_64::_mm256_set_m128i(masks[0].0, masks[1].0);
            let masked = x86_64::_mm256_and_si256(self256, mask256);

            let eq = x86_64::_mm256_cmpeq_epi8(masked, x86_64::_mm256_setzero_si256());
            let mv = x86_64::_mm256_movemask_epi8(eq) as u32;
            let e0 = MASKED_VALUES[15 - (mv ^ 0xffff_ffff | 0x0001_0000).leading_zeros() as usize];
            let e1 = MASKED_VALUES[31 - (mv & 0xffff ^ 0xffff | 0x0001).leading_zeros() as usize];

            let m = masked;
            let m = x86_64::_mm256_or_si256(m, x86_64::_mm256_srli_epi16::<1>(m));
            let m = x86_64::_mm256_or_si256(m, x86_64::_mm256_srli_epi16::<2>(m));
            let m = x86_64::_mm256_or_si256(m, x86_64::_mm256_srli_epi16::<4>(m));
            let m = x86_64::_mm256_or_si256(
                x86_64::_mm256_srli_epi16::<1>(m),
                x86_64::_mm256_set_epi64x(e0.1, e0.0, e1.1, e1.0),
            );
            let ret = x86_64::_mm256_andnot_si256(m, mask256);
            Self(x86_64::_mm_or_si128(
                x86_64::_mm256_castsi256_si128(ret),
                x86_64::_mm256_extracti128_si256::<1>(ret),
            ))
        }
    }
}

これでとりあえずは同時に2方向ずつの飛び利きを求めることはできた。

AArch64

x86_64用のSIMDを頑張って実装したものの、じつは手元のメインの開発機としては最近M1 MBPを使っている。

のでこちらでもSIMDを使えるようにしたいと思い、NEONを使って同様に128bit計算を同時に行うように実装してみた。

use std::arch::aarch64;

pub(crate) struct Bitboard(aarch64::uint64x2_t);

SIMD Bitboard for AArch64 NEON by sugyan · Pull Request #19 · sugyan/yasai · GitHub

だいたいx86_64と同じように出来るかな…と思ったが まぁまぁ違うところがあって苦戦した。

同値判定、ゼロ値判定

x86_64::_mm_test_all_zerosi32 で値を返してくれたが、 aarch64 にはそれに相当するようなものは無い。 aarch64::vceqq_u64 なり aarch64::veorq_u64 なりで2値の比較をしても、いちいちその結果から中身を取り出して値を確認しないと bool を返すことはできない。

検索して調べてみたが、どうやら aarch64::vqmovn_u64 を使って3命令でゼロ値判定するのが最良のようだった。

impl Bitboard {

    ...

    pub fn is_empty(&self) -> bool {
        unsafe {
            let vqmovn = aarch64::vqmovn_u64(self.0);
            let result = aarch64::vreinterpret_u64_u32(vqmovn);
            aarch64::vget_lane_u64::<0>(result) == 0
        }
    }
}

VQMOVN (Vector Saturating Move and Narrow) copies each element of the operand vector to the corresponding element of the destination vector. The result element is half the width of the operand element, and values are saturated to the result width.

https://developer.arm.com/documentation/dui0489/h/neon-and-vfp-programming/vmovl--v-q-movn--vqmovun

とのことで、bit幅を半分にし 値がオーバーフローする場合は飽和させてコピーする、というものらしい。 以下の記事が詳しくて分かりやすかった。

qiita.com

ゼロ値との比較に使うぶんには飽和の影響は無く、正しく 0 か否かを保持したまま uint64x2_t から uint32x2_t に変換できる、ということになる。あとはこれを uint64x1_t として解釈した上でその要素として単一の u64 を取り出して 0 と比較すれば良い。

飛び利き計算

香車の利きはx86_64と同様にQugiyの手法で問題なかったが、飛車・角行の利きには同じことをしようにも x86_64::_mm256_movemask_epi8 と同等なものが見つからなかった。 各8bitから値をかき集めて単一の u32(もしくは u16)を得る、のようなことが簡単には出来ないらしい。

また、 x86_64::_mm_slli_si128 のようにlaneを跨いだbit shiftも出来ないようだった。まぁdecrementの桁借り計算で必要なのは64bitのlane丸ごと移動だから aarch64::vextq_u64 でどうにかはなるのだけど…

ここはまぁ仕方ないか、ということで前回の記事で書いた Leading/Trailing Zeros を使う手法で実装した。毎回 [u64; 2] の値を取得して分岐して…というのが必要になるのは微妙ではあるが そこまで遅いというわけでもないと思う。

use std::mem::MaybeUninit;

const MASKED_VALUES: [[u64; 2]; Square::NUM + 2] = {
    let mut values = [[0; 2]; Square::NUM + 2];
    let mut i = 0;
    while i < Square::NUM + 2 {
        let u = (1_u128 << i) - 1;
        values[i] = [u as u64, (u >> 64) as u64];
        i += 1;
    }
    values
};

impl Bitboard {

    ...

    #[inline(always)]
    fn values(self) -> [u64; 2] {
        unsafe {
            let m = MaybeUninit::<[u64; 2]>::uninit();
            aarch64::vst1q_u64(m.as_ptr() as *mut _, self.0);
            m.assume_init()
        }
    }
    fn sliding_positive(&self, mask: &Bitboard) -> Bitboard {
        let m = (*self & mask).values();
        let tz = if m[0] == 0 {
            (m[1] | 0x0002_0000).trailing_zeros() + 64
        } else {
            m[0].trailing_zeros()
        };
        Self(unsafe {
            aarch64::vandq_u64(
                mask.0,
                aarch64::vld1q_u64(MASKED_VALUES[tz as usize + 1].as_ptr()),
            )
        })
    }
    fn sliding_negative(&self, mask: &Bitboard) -> Bitboard {
        let m = (*self & mask).values();
        let lz = if m[1] == 0 {
            (m[0] | 1).leading_zeros() + 64
        } else {
            m[1].leading_zeros()
        };
        Self(unsafe {
            aarch64::vbicq_u64(
                mask.0,
                aarch64::vld1q_u64(MASKED_VALUES[127 - lz as usize].as_ptr()),
            )
        })
    }
}

Iterator

しかしどうも実装して実際にbenchmarkを取ってみると、すごく遅い。SIMD使わない方が速いくらい。一体何故…?と思って調べてみると、Bitboardの計算自体ではなく、そこから Square を取り出す処理が遅くなっているようだった。

合法手の生成にはそれぞれ利きの通る位置への移動を列挙したりするので、論理演算から得たBitboardからすべての対応する Square を列挙するような操作が多くなる。 これには Iterator を実装し その内部で fn pop(&mut self) -> Option<Square> を呼ぶのが初期の実装だった。

impl Bitboard {
    pub fn pop(&mut self) -> Option<Square> {
        let mut m = unsafe {
            let mut u = std::mem::MaybeUninit::<[u64; 2]>::uninit();
            aarch64::vst1q_u64(u.as_mut_ptr() as *mut _, self.0);
            u.assume_init()
        };
        if m[0] != 0 {
            unsafe {
                let sq = Some(Square::from_u8_unchecked(Self::pop_lsb(&mut m[0]) + 1));
                self.0 = aarch64::vsetq_lane_u64::<0>(m[0], self.0);
                sq
            }
        } else if m[1] != 0 {
            unsafe {
                let sq = Some(Square::from_u8_unchecked(Self::pop_lsb(&mut m[1]) + 64));
                self.0 = aarch64::vsetq_lane_u64::<1>(m[1], self.0);
                sq
            }
        } else {
            None
        }
    }
    fn pop_lsb(n: &mut u64) -> u8 {
        let ret = n.trailing_zeros() as u8;
        *n = *n & (*n - 1);
        ret
    }
}

impl Iterator for Bitboard {
    type Item = Square;

    fn next(&mut self) -> Option<Self::Item> {
        self.pop()
    }
}

内部の値の trailing_zeros() で rightmost set bit の位置を探して Square を作り、その後そのbitだけを 0 にclearする。これを値がすべてゼロになるまで繰り返す。

とはいえ trailing_zeros() を得るのはSIMDでは出来ないのでどうしても一度 u64 の値で読み込む必要はある。各要素で分岐は必要だし、それを pop_lsb で bit clear した値を書き戻す必要もある。

どうもこの読み込みと書き込みを繰り返す操作が、x86_64のときはそこまでオーバーヘッドを気にすることが無かったが AArch64では特に遅くなっていて影響が出ているようだ。 合法手生成でいうと持駒を打つ場所を列挙するときなど、多くの候補がある(すなわち多くの 1 が存在している)Bitboardに対する Square 列挙の繰り返しで影響が大きくなっていることが分かった。

では、SIMDレジスタからの読み込みを1回だけにして、その後はSIMDレジスタを使わずに処理するようにしたらどうか?ということで Iterator 用のstructを用意し、 IntoIterator を実装してそちらに pop の処理を任せる実装にしてみた。

pub(crate) struct SquareIterator([u64; 2]);

impl SquareIterator {
    #[inline(always)]
    fn pop_lsb(n: &mut u64) -> u8 {
        let pos = n.trailing_zeros() as u8;
        *n &= n.wrapping_sub(1);
        pos
    }
}

impl Iterator for SquareIterator {
     type Item = Square;

     // #[inline(always)]
     fn next(&mut self) -> Option<Self::Item> {
        if self.0[0] != 0 {
            return Some(unsafe { Square::from_u8_unchecked(Self::pop_lsb(&mut self.0[0]) + 1) });
        }
        if self.0[1] != 0 {
            return Some(unsafe { Square::from_u8_unchecked(Self::pop_lsb(&mut self.0[1]) + 64) });
        }
        None
    }
}

impl IntoIterator for Bitboard {
    type Item = Square;
    type IntoIter = SquareIterator;

    fn into_iter(self) -> Self::IntoIter {
        unsafe {
            let m = std::mem::MaybeUninit::<[u64; 2]>::uninit();
            aarch64::vst1q_u64(m.as_ptr() as *mut _, self.0);
            SquareIterator(m.assume_init())
        }
    }
}

これで、以前と同様に for sq in bb { ... } のようにfor loopですべての Square 列挙の繰り返しができるし、SIMDレジスタからの読み込みは1回で済むから速くなるはず… と思ったがbenchmarkとってみると効果が無い。。

ところが、上述のコードで // #[inline(always)]コメントアウトしている部分、fn next(&mut self) -> Option<Self::Item>inline(always) attribute をつけると劇的に改善された。どうして…!?

一応、手元の環境で様々な実装を比較しながらbenchmarkを回してみた。

https://gist.github.com/sugyan/03c1e73a253b4591dc9f29f9609ad93e

SIMDを使わずに [u64; 2]Iterator専用structを用意するのが、やはり素朴で最も速いようだ。それと同様の実装をしているはずの aarch64_iterator, x86_64_iterator や、他の pop を呼ぶ実装などはすべて遅い。 同じ実装だが next(&mut self)inline(always) をつけただけの差異しか無い aarch64_iterator_inline, x86_64_iterator_inline は最速に近い値になった。これだけでこんなに変わるものなのか…

https://godbolt.org/z/r3qj898rE

inline(always) の有無による違いも見てみたけどイマイチよく分からない… どういうことなんだろう。

ともかく、これをつけるかどうかで3〜4倍も速度差があることが分かった。x86_64では2.5倍くらいの差しか出ていなくて、この影響に気付けなかったようだ。 とはいえIterator専用structを用意する手法の方が良さそうなのでどちらもこれを採用することにした。

WebAssembly

ここまでやったらwasm32の 128-bit packed SIMD にも対応してやるぞ、ということで同様にやってみた。

use std::arch::wasm32;

pub(crate) struct Bitboard(wasm32::v128);

Implement wasm32 simd128 Bitboard by sugyan · Pull Request #23 · sugyan/yasai · GitHub

std::arch::wasm32APIが扱いやすく、それほど苦労することなく実装できた。 同値判定、ゼロ値判定には wasm32::u64x2_all_truewasm32::v128_any_true が使えるので便利。 SIMDレジスタ間の読み書きオーバーヘッドなどは特に考える必要なさそうで、とにかく命令数を少なめにしておけば良い感じ。素朴に読み込んで Leading/Trailing Zeros を使って飛び利きを計算するようにした。

各moduleで unit test を書いてあるがwasm32の場合の実行の仕方がよく分からなかった。とりあえず Wasmer を入れて、 wasm32-wasi でbuildして実行するようにしてみた。

export RUSTFLAGS="-C target-feature=+simd128" 
cargo clean
cargo build --tests --target=wasm32-wasi --features simd
wasmer run target/wasm32-wasi/debug/deps/yasai-*.wasm

また、これを使って wasm-packSIMDを有効にしたアプリケーションを作ろうとすると wasm-opt が古くてSIMDに対応していないらしく、自前で最新版をinstallして使う必要があるようだった。

ハマりどころはそれくらいだっただろうか。

Benchmark

実際どれくらい効果が出るか、と手元のPCで perft 5 を計測してみると。

x86_64

MacBook Pro (13-inch, 2017, Four Thunderbolt 3 Ports) 3.3 GHz Dual-Core Intel Core i5

$ cargo +nightly bench perft::bench_perft_5

...

test perft::bench_perft_5_from_default       ... bench: 357,959,818 ns/iter (+/- 5,399,589)

$ cargo +nightly bench --features simd perft::bench_perft_5

...

test perft::bench_perft_5_from_default       ... bench: 276,508,910 ns/iter (+/- 1,445,515)

NPS約1.3倍。

AArch64

MacBook Pro (14-inch, 2021) Apple M1 Pro

$ cargo +nightly bench perft::bench_perft_5

...

test perft::bench_perft_5_from_default       ... bench: 194,382,054 ns/iter (+/- 2,067,121)

$ cargo +nightly bench --features simd perft::bench_perft_5

...

test perft::bench_perft_5_from_default       ... bench: 163,265,674 ns/iter (+/- 4,929,186)

NPS約1.2倍。

WebAssembly

これはWebで確かめられるので実際にSafari以外のブラウザで確かめてみていただきたい。 (SafariSIMD対応していなかったの知らなかった…)

手元のChromeでは1.5倍ちかく速くなる。

感想

Benchmarkの数値を見てわかる通り、古いPCでSIMDつけてちょっと速くするより 良いPCに替えた方がSIMD無しでも圧倒的に速いので、まぁ結局はマシンパワーが強い方がつよい。

とはいえこういったBitboardの実装を差し替えるだけで効果が出る高速化というのはやる価値はちゃんとあると思うのでやってみて良かった。勉強になることもたくさんあった。

逆にこのBitboard関連ではこれ以上の高速化はそんなに望めないので、合法手生成をもっと速くしようと思ったらあとは生成ロジックの方をちゃんと見直す必要がある、ということになる。今後はここを頑張っていきたい。

Bitboardでleading/trailing zerosを使って飛び駒の利きを求める

前置き

Rustで将棋の高速合法手生成ライブラリを作り始めていて、

その後 shogi_core というライブラリが出来てきたのでそれを使うように変更したりした。

それまでは基本的に Rust版Apery の実装を参考にしていたが、駒の利きを求めるところで改善したいところがあり、Qugiyの手法を取り入れてみようと検討した。

そこで色々試したり考えているうちに別のやり方を思いついたのでメモしておく。

先にオチを書いておくと、チェスの世界では既出の手法でした。

飛び利き計算いろいろ

基礎知識

基本的にBitboardを使っての駒の利きは表引きすることで求める。 事前に「この位置にいるこの手番のこの駒は、ここに動ける」というのを全パターン計算してあれば、テーブルから引いてくることで一発で移動先が求まる。

struct PieceAttackTable([[Bitboard; Color::NUM]; Square::NUM]);

impl PieceAttackTable {
    fn new() -> Self {
        ... // 事前計算
    }
    pub fn attack(&self, sq: Square, c: Color) -> Bitboard {
        self.0[sq.array_index()][c.array_index()]
    }
}

これを全駒種に対して用意しておけば良い。 …のだが、香車・角行・飛車の3種の「飛び駒」(という名称が正しいのかよく分かっていない) についてはそうはいかない。これらは「最初に他の駒にぶつかるまで一方向にどこまでも動ける」ので、利きの届く範囲は他の駒の位置に依存する。 ので簡単な表引きでは求めることができず、様々な手法が考案されてきていた。

以下、利きを求めるために使われる「盤上に(自陣敵陣問わず)駒が存在している位置」を表すBitboardをOccupancyと呼ぶことにする。

縦型Bitboardでの香車の利き

利きはOccupancy依存とは言っても、考慮すべきは進む方向のSquareだけで良い。一段目にいる後手香車が駒の存在を意識すべきは同じ筋の二段目〜八段目だけ(九段目は居ても居なくてもそこまで届く、と考えられる)なので、7つのSquareにそれぞれ駒が居るか居ないかの情報だけあれば良い。 ということはその7つのSquareの状態は全パターン列挙しても高々128(1 << 7)通りなので、それくらいならテーブルを用意しておいて引ける、という考え方。

pub struct LanceAttackTable(Vec<Vec<Vec<Bitboard>>>);

impl LanceAttackTable {
    const MASK: usize = (1 << 7) - 1;

    fn new() -> Self {
        let mut table = vec![[[Bitboard::empty(); 1 << 7]; Color::NUM]; Square::NUM];
        ... // 事前計算
    }
    pub fn attack(&self, sq: Square, c: Color, occ: &Bitboard) -> Bitboard {
        let index = Self::MASK & ((occ.value() >> ((sq.file() * 9) + 1)) as usize);
        self.0[sq.array_index()][c.array_index()][index]
    }
}

少しテーブルのサイズは大きくなるが、それほど影響ない範囲で事前計算しておける。

ここで、縦型にレイアウトしているBitboardだとOccupancyの状態からindexに変換する際に利点が出てくる。 縦の動きの場合は対象の筋のbitが連続して並ぶので、これを適当にshiftしてmaskかければパターンのindexを簡単に求めることができる。

角行・飛車の利き (Magic Bitboard / PEXT)

そして角行・飛車でも同じように全パターンを列挙して表引きしたいが、香車の場合と違って

  • パターンの数が多く巨大なテーブルが必要、事前計算にも時間がかかる
    • 具体的には 角行が20,224パターン、飛車が495,616パターン
  • 対象bitが連続して並ばないのでindexの計算が容易ではない

といった問題が出てくる。 前者については apery_rustonce_cell crateを使うことで最初に呼ばれたときにだけ初期化の計算する方式をとっている。

後者について。 前述の通り、一直線に並んでいる香車のような動きなら全体をshiftすることで簡単にindexを計算できるが、横や斜めの動きに関しては注目すべきbitが飛び飛びになるので同じようにはできない。 何とかして対象となるmaskの位置のbitだけの値をかき集めてきて数値を得たい。

 .##.#.#. ##.#.#.# | occupancy
 ..#..#.. #..#..#. | mask
(  1  0   1  1  0 )
----------------------
             10110 | f(occupancy, mask)

これをまさしく実現してくれるのが BMI2のPEXT命令。これを使うことで香車の利きのindex算出と同様のことができるので、一発で表引きできるようになる。

ただしPEXTは64bit整数にしか適用できないので、Bitboard内部で2つの64bit整数に分けて状態保持している場合は p[0] | p[1] のようにmergeした上でPEXTにかけるなどの工夫が必要になる。この方法が問題なく動作するためには縦型でも1〜7筋と8,9筋で分ける必要が出てくる、など制限が生まれる。

歴史的にはこのPEXT方式が出てくる前にMagic Bitboardが発明されていた、のかな。PEXTと同様の処理をMagic Numberを乗算することで実現する、というもの。

PEXTの方がMagic Numberを用意する必要がないという利点がある一方、

  • 特定の命令セットに依存することになる
    • 全環境での動作をサポートするには何らかのfallbackは必要
  • 一部のCPUアーキテクチャでは異常に遅いらしい

といった問題はある。また巨大なテーブルが必要という問題は依然として残り続ける。 このあたりは以下の記事などのBitboard関連シリーズがとても詳しく、勉強になりました。

yaneuraou.yaneu.com

Qugiyの手法

近年の将棋ソフトウェアではおそらく前述のMagic Bitboard/PEXTでの表引きが主流となっている(要出典)が、2021年のWCSC31で登場したQugiyが新しい手法を考案して話題になったらしい。以下の記事が詳しい。

yaneuraou.yaneu.com

要するに加減算での桁上げ/桁借りの伝搬を利用して最初に立っているbitまでのmaskを作る、というのが基本概念のようだ。

 (.##.#.#. ##.#.... | occupancy)
  #..#.#.# ..#.#### | !occupancy
  #..#.#.# ..##.... | !occupancy + 1
---------------------  
  ............##### | !occupancy ^ (!occupancy + 1)

減算ならより簡潔で

  .##.#.#. ##.#.... | occupancy
  .##.#.#. ##..#### | occupancy - 1
---------------------  
  ........ ...##### | occupancy ^ (occupancy - 1)

この概念でいうと飛び飛びのmaskになる角行・飛車でも同様に考えることができ、

  .##.#.#. ##.#.... | occupancy
  ..#..#.. #..#..#. | mask
---------------------  
  ..#..... #..#.... | masked = occupancy & mask
  ..#..... #...#### | masked - 1
---------------------
  ........ ...##### | masked ^ (masked - 1)
  ..#..#.. #..#..#. | mask
---------------------
  ........ ...#..#. | mask & (masked ^ (masked - 1))

と、mask上で最初にoccupancyにぶつかるまでの利きをbit演算だけで求めることができる。

ただし、この方法は桁上げ/桁借りが伝搬する方向にしか適用できず、つまり上から下にレイアウトされている縦型Bitboardでは下方向か左方向にしか使えない。ので、上方向や右方向に対してはbyte swapやunpackを使うなどの工夫をしている、とのこと。

あとは4方向にそれぞれ別々に求める必要があるところを、左下方向の2つと右上方向の2つそれぞれをAVX2命令で並列処理するようにしているらしい。

Leading/Trailing Zerosを使う

Qugiyの手法を実装しようとしたものの、もうちょっと簡単に同じようなことを実現できないかな、と考えていて、この方法を思いついた。

「要するに加減算での桁上げ/桁借りの伝搬を利用して最初に立っているbitまでのmaskを作る」というのが肝だったので、つまり最初に立っているbit(rightmost set bit)を見つければそこまでのmaskを作ってかけてやれば良いだけではないか、と。

  .##.#.#. ##.#.... | occupancy
  ..#..#.. #..#..#. | mask
---------------------  
  ..#..... #..#.... | masked = occupancy & mask
  ........ ...#.... | r = rightmost set bit of masked
  ........ ..#..... | r << 1
  ........ ...##### | (r << 1) - 1
  ..#..#.. #..#..#. | mask
---------------------
  ........ ...#..#. | mask & ((r << 1) - 1)

結局、知りたいのはoccupancyなりoccupancy & maskなりの「最初に立っているbitの位置」だった。これは最下位から見る場合はBitboardからSquareを取り出すときにも使われている Trailing Zeros そのもの。

で、同様のことは逆方向でも考えることができて、こちらは Leading Zeros を使ってleftmost set bitを求めることが出来る。作るmaskが反転する必要あるだけ。

  .#..###. ##.#.... | occupancy
  ..#..#.. #..#..#. | mask
---------------------  
  .....#.. #..#.... | masked = occupancy & mask
  .....#.. ........ | l = leftmost set bit of masked
  ......## ######## | l - 1
  ######.. ........ | !(l - 1)
  ..#..#.. #..#..#. | mask
---------------------
  ..#..#.. ........ | mask & !(l - 1)

ということで、実はこんなコードで簡単に両方向の利きを計算することが出来た。

const BB_1A: Bitboard = Bitboard::single(Square::SQ_1A);
const BB_9I: Bitboard = Bitboard::single(Square::SQ_9I);

fn sliding_positive(occupancy: Bitboard, mask: Bitboard) -> Bitboard {
    let tz = (occupancy & mask | BB_9I).to_u128().trailing_zeros();
    mask & unsafe { Bitboard::from_u128_unchecked((1 << (tz + 1)) - 1) }
}
fn sliding_negative(occupancy: Bitboard, mask: Bitboard) -> Bitboard {
    let lz = (occupancy & mask | BB_1A).to_u128().leading_zeros();
    mask & unsafe { Bitboard::from_u128_unchecked(!((1 << (127 - lz)) - 1)) }
}

occupancy & mask がemptyになるとleading/trailing zeroの値がはみ出てしまうので便宜上 端っこのbitを立ててしまって、それぞれ正方向なら .trailing_zeros()、負方向なら .leading_zeros() で最初に立っているbit位置を求められる。あとはその位置から両端までのmaskを作って元のmaskと論理積とったものを返すだけ。両端までのmaskは事前に作っておいて表引きするだけでも良いかもしれない。

ソフトウェア的に実装すると .leading_zeros().trailing_zeros() と比較してコストがかかる処理になったりもするらしいが、プロセッサよっては lzcnt など専用命令があったりもするのでコンパイルオプションなどを適切に使用すればそれほど気にしなくても良さそう…?

実験と考察

ということで、巨大テーブルも必要なく、PEXTのような特定の命令セットに依存することもなく、それほど負荷がかからない処理で飛び利きを計算することが出来た。

自分のライブラリで実装してみたところ、PEXTでの表引きと比較すると少し遅くなる、程度の差異となった。

github.com

やはり4方向それぞれ計算することになるのはPEXTの表引き一発と比較すると遅くはなってしまうと思われる。巨大テーブルが不要になり多少時間かかっていた初期化処理をする必要が無くなったのは恩恵だが…

Qugiyと同様に、同じ方向のものはSIMDで同時に計算する、などの工夫を入れることで多少は改善されるだろうか。

ただこの例ではshogi_coreのインタフェースとして内部の値に触れるのがu128だったからこうしているが、SIMDを使って実装しているとleading/trailing zerosは簡単には求められないかもしれない…。そのへんの相性の悪さもあって今まで採用されてこなかった、とかもあるのだろうか…?

この手法を思いついて、今までこういう手法なかったでしょうか?と開発者Discordで尋ねたところ最初に書いた通りチェスの世界では "Classical Approach" として知られていたものだった。 実装が簡単な反面、速度的には限界があるのかもしれない。

とはいえそこまですごい遅いわけでもなく簡単に実装できる手法ではあるので、知っておいて損は無いと思う。Classicalではあったが自分でこの手法を思いつくことができたのは良かった。