60bit整数を高速に素因数分解する

Library CheckerFactorizeという問題でfastestを取ったので、そのコードの解説です(速度ランキング)。 Factorizeは、 10^{18}までの整数が高々100個与えられ、それらを素因数分解せよ、という問題です。

概要

上位陣(13ms~15ms)の解法では、ECM楕円曲線を用いた素因数分解法)が使われているようですが、よくわからないのでこれには手を出しません。 ECMを用いていない提出の中で最速の提出はhqztrue氏の16msであり、20%以上の差があります。 この提出では、ポラードのロー法という手軽に実装できるアルゴリズムが用いられていました。 これを効率化してみます。

ポラードのロー法には、乱数を連鎖的に生成する部分が存在します。 この部分がレイテンシに支配されていて並列性に乏しいという問題点に着目しました。 これに対し、異なる二つの系列の乱数列を同時に生成する手法を導入します。 これにより、プログラムの並列度をスループットに支配されるところまで上げることができます。 並列化のオーバーヘッドがかなりあり、二倍速くなるというわけではないですが、それでも12msと上位陣の解法を上回る速度を達成できました。

基礎知識

純粋にアルゴリズムの基礎知識です。 CPUの特性を生かした高速化に関する記述はほぼないので、高速化テクニックが知りたい場合は「本題:私の提出コードの解説」まで読み飛ばしても多分大丈夫です。 ポラードのロー法を“貼って”いるけれど、何をやっているのかあまりよく知らない、という場合は読んでみると面白いかもしれません。

ポラードのロー法

ポラードのロー法は、合成数 nを入力すると、その非自明な因子( n 1 nでない約数) aを返す乱択アルゴリズム(ラスベガス法)です。 ここで、返ってくる a素数とは限らない点に注意してください。 とはいえ、素数判定と組み合わせ、素数でなければ再帰的にロー法を適用する、とすることで素因数分解を行うことができます。

乱数を使って素因数分解するというのはどういうことでしょうか。 単純な方法として、 [0, n)の範囲の整数を一様ランダムに生成し、それと nの最大公約数(greatest common divisor, GCD)を求める、という方法があります。 これを繰り返していれば、そのうち aが返ってくることがあるでしょう。 しかし、これには平均で aに比例した回数の乱数生成が必要です。

誕生日のパラドックスを考えると、次のような方法で乱数生成回数を減らせます。  [0, n)の範囲の整数をランダムに \sqrt a個程度生成し、それら同士の差を計算します。 ここで計算された「差」一つ一つに対して、 nとの間の最大公約数を求めれば、乱数生成は \sqrt aに比例した回数で済みます。 しかし、乱数を二つ取り出す組み合わせの数を考えると、「差」は a/2個ほどあるので、計算量はやはり aに比例してしまいます。 なんとか無意味な組み合わせの「差」を求めることを回避できないでしょうか。

これに対してはうまい方法があり、それがポラードのロー法です。 ここまでは一様ランダムな真の乱数を考えていましたが、そうではなく、あえて疑似乱数を使うというものです。 例えば、なんらかの整数 Cに依存した疑似乱数系列 x_i \mapsto x_i^2 + C \bmod nがよく用いられます。‘

「差」と nの間の最大公約数を求めたときに aが返ってくるのは、「差」を取った二つの乱数が \bmod aで一致していた場合なので、乱数の aでの剰余に着目します。 いきなり aでの剰余に着目してもいいのですが、それだと当たり前のことを言っているだけのようにも感じてしまうかもしれません。 そこで、あえて任意の Kでの剰余の場合にどうなるかを考え、その後 K = aという特殊な状況ではどのように性質が変わるかを見てみます。

さて、この疑似乱数系列を \mod\! Kで見た y_i = x_i \bmod Kを考えます。 この時、 y_i \mapsto y_{i+1}という変換は一意に決まる(写像になる)でしょうか?  Kがなんでもない数であれば、これは一意に決まるとは限りません。  x_i y_iよりも多くの情報を持っているので、 y_iが決まっても y_{i+1}には複数の可能性があります。

例で考えてみます。  x_i \mapsto x_i^2 + 1\bmod 15という疑似乱数生成式は、  x_1 = 1に対して1→2→5→11→2→5→11→……という疑似乱数系列を生み出します。 ここで、 K = 6で割った余りを考えると、1→2→5→5→2→5→5→……となっています。 つまり、5の次が5だったり2だったりしています。  y_i \mapsto y_{i+1}という変換が一意でないということです。

しかし、 K nの約数の場合は状況が違います。 たしかに x_i y_iより多くの情報を持っていますが、それは y_i y_{i+1}の関係には影響を及ぼしません。 なぜなら、 y_i \mapsto y_i^2 + C \bmod Kが成り立つからです( \bmod rsで成り立つ式は \bmod rでも成り立つ、という性質を使いました)。 よって、  K nの約数の場合は変換が一意です。 変換が一意なので、系列中に同じ値が出現すると、以降は全く同じ系列が繰り返されることになります(上記例では5→5となったからといって、その後も5→5→5→5→……となるわけではないのと対比してみてください)。 以下では、系列が繰り返される周期を Pと表記します。

話を戻すと、乱数 x_iと乱数 x_jの差 x_i - x_jを求め nとの最大公約数を求める手順で、無意味な組み合わせを回避したいという話でした。 無意味な組み合わせとは、 x_i - x_j aの倍数にならないことが明らかな組み合わせのことです。 これは、 K=aで考えれば、 y_i y_jが必ず異なることが明らかな組み合わせのことです。 ポラードのロー法では、疑似乱数系列に周期があることをうまく活用し、無駄な組み合わせを回避します。 例えば、 y_{i+1} y_{i+3}が異なっていれば、 y_i y_{i+2}も異なっています(対偶を考えれば明らかです)。 よって、 x_{i+1}-x_{i+3} nの最大公約数が1だった場合、 x_i-x_{i+2} nの最大公約数を求める必要はありません。

この考え方、つまり添え字の差が同じペアは調べても無駄、という事実を活用すると、調べるべきペアの個数を P個程度に減らすことができます。  x_i x_{i+1}の差と nの最大公約数を求める、  x_i x_{i+2}の差、 x_i x_{i+3}の差、 x_i x_{i+4}の差、……とP回繰り返せばよいからです(※)。 乱数を P個生成した時、ペアは P^2/2個程度ありますが、調べる必要があるのは P個程度でしかない、ということです。

周期 P K=a以下*1であることが明らかですが、 y_{i+1} = x_{i+1} \bmod K = y_i^2 + C \bmod Kという疑似乱数生成式を使った場合、周期 P \sqrt K=\sqrt a程度になることが経験的に知られています。  nの任意の因子 cについて、 \sqrt c個程度のペアを調べると cが見つかるということなので、基本的には小さい因子から見つかっていきます。  nの最大の素因子はこの方法で見つけなくていい( nを見つかった因子で割るということを繰り返して残った数が素数かを判定すればよい)ので、 nの二番目に大きな素因子の大きさの正の平方根に比例するくらいの時間がかかるとみてよいでしょう。  nの二番目に大きな素因子の大きさは高々 \sqrt nなので、 n^\frac14に比例した回数の乱数生成により素因数分解できることが経験的に知られています。

※実際には、 iをどう選ぶかという問題があります。  y_1が循環に含まれるとは限らないからです(例えば、上で紹介した1→2→5→11→2→5→11→……の例では、1は循環に含まれていません)。 これは、適宜 i自体を進めていく、という方法で解決可能です。 ポラードが提案した進め方(二回に一回進める)よりも乱数生成回数が少なくなる進め方( 2^t回に一回、系列を 2^{t-1}だけ進める)をブレントが提案しています。

参考

以下の記事に、アニメーションを用いたわかりやすい解説があります。

ポラード・ロー素因数分解法について #競技プログラミング - Qiita

モンゴメリ乗算

 \mod\!nの計算は、普通にやると除算命令が必要です。 除算命令は一般に遅いのでどうにか回避したいです。

モンゴメリ乗算は、バレット還元と同様に、いくらかの事前計算を行うことで、剰余演算を除算命令無しに行う方法です。 モンゴメリ乗算では、モンゴメリ表現と呼ばれる数値表現で値を保持したまま計算が進みます。 モンゴメリ表現は、 \bmod nの代表元 Xの代わりに、 2^{64}などの切りのいい数を掛けた 2^{64}X \bmod nを使います。

事前計算で求めるものは、モンゴメリ表現から通常の表現に戻すときに使う値 Rです(他にも用途ごとに他の数を事前計算することがあります)。 これは \bmod 2^{64}での nでの逆元であり、  nR \equiv 1 \mod 2^{64}を解けば求まります。 この Rを用いると、乗算を以下のように行うことができます。

元の表現 モンゴメリ表現
 X  A = 2^{64}X \bmod n
 Y  B = 2^{64}Y \bmod n
 XY  C = 2^{64}XY\bmod n
 \;\;\;\:= 2^{-64}AB\bmod n

さて、この 2^{-64}AB\bmod nはどのように求めたらよいでしょうか? ここで \bmod nが発生していたら意味がありません。 これを事前計算された Rを使って回避するのが、モンゴメリ乗算です。  Z \equiv 2^{-64}D \mod nとなる Zは次のように計算できます。

 Z = \frac{D}{2^{64}} - \frac{(DR \bmod 2^{64})n}{2^{64}}

この手順には 2^{64}での除算や 2^{64}での剰余を求める部分がありますが、64bitCPUであれば実質ノーコストである点に注意します。

なぜこれでいいのでしょうか?まず、 Zが整数であることから確認します。

 D - (DR \bmod 2^{64})n \equiv D - DRn \equiv D - D\times1 \equiv 0 \mod 2^{64}

なので、 D - (DR \bmod 2^{64})n 2^{64}の倍数であり、 2^{64}で割っても整数です。 また、  D (DR \bmod 2^{64})n 2^{64}の倍数とは限りませんが、その差が 2^{64}の倍数であるということは、下位64bitは同じになっています。 よって、もし \frac{D}{2^{64}} \frac{(DR \bmod 2^{64})n}{2^{64}}を切り捨て除算で計算しても、結果は同じになります。  Z = \frac{D - (DR \bmod 2^{64})n}{2^{64}}のように計算すると D - (DR \bmod 2^{64})nの下位64bitの計算が一見必要そう(下位64bitは除算で消えるにせよ、ボローが生じるかを計算しないといけなさそう)ですが、実はそれは不要、ということを意味しています。

つぎに、 Z \equiv 2^{-64}D\mod nであることを確認します。

 Z = 2^{-64}D - 2^{-64}(DR \bmod 2^{64})n \equiv 2^{-64}D \mod n

大丈夫そうです。

このように、 2^{64}などの切りのいい数を掛けた表現で持つことで、多少の事前計算と引き換えに剰余演算を除去するのが、モンゴメリ乗算です。

なお、単に 2^{64}を掛けた表現で持っているだけなので、加減算は普通の剰余環でのものと同じようにして実行可能です。

ポラードのロー法では非常に多くの剰余演算が使われるため、高速化のためには事実上必須のテクニックです。 以下では、ポラードのロー法といったとき、暗黙でこのテクニックが適用されているものとします。

素数判定

ポラードのロー法で素因数分解するためには、素数であるかを(素因数分解によらずに)判定する必要があります。 いろいろな高速な手法が提案されていますが、素因数分解にかかる時間が圧倒的なため、素数判定の高速化に躍起になる必要はありません。 今回はミラーラビン素数判定法を愚直に実装したものを使いました。 ミラーラビン素数判定法は確率的判定法ですが、 2^{64}までの整数の素数判定に対しては、特定の7つの底を使えば間違いなく判定できることが知られています。

Binary GCD

Binary GCD(Stein's algorithm)は、最大公約数を求めるアルゴリズムです。 二進法を採用している計算機(今日のほぼすべての計算機)に適しており、高速に動作します。

最大公約数を求めるアルゴリズムとして有名なユークリッドの互除法は、剰余演算が含まれるため、あまり高速ではありません。 アルゴリズムの進行とともに入力依存で法が変化するため、剰余演算を事前計算で高速化することも不可能です。 Binary GCDは除算を含まないため、高速に動作します。

以前の記事(最大公約数をもっと高速に求める番外編(その2) - よーる)では、C言語コードを工夫して書けばclangが最適なコードを吐くようにできることを紹介しました(インラインアセンブリに手を出さなくて済むということです)。 当時はgccだと最適なコードとならない問題がありましたが、現在のgcc(ジャッジサーバーで使われているのはgcc 13.2)ではほぼ最適なコードが出るようになっているようです。

60bit整数同士の最大公約数を求める場合を実験してみると、このコードのループは平均で40回程度実行されるようです。 Binary GCD (GNU MP 6.1.2)によれば1bitあたり0.68回の反復が必要、と書かれているのでこれと合致します。 ループ一周は5cycleで実行可能*2なため、最大公約数を求めるのにかかる時間は200cycle程度です。

ポラードのロー法で、最大公約数をまとめて計算するテクニック

ポラードのロー法を愚直に実装すると、「乱数を生成し、差を計算し、それと nの最大公約数を求める」というのを繰り返すコードになると思います。 しかし、これはかなり遅くなります。 先ほど示したように、最大公約数を求めるのにかかる時間は200cycle程度です。 これは乱数の生成にかかる時間が10cycleや20cycleであるのに対して10倍も大きいです。 これではアルゴリズムの実行時間の大半が最大公約数を求めることに費やされてしまいます。

これには面白い解決策があります。  nとの間の最大公約数を求める数を、1つの差とするのではなく、 M個の差を全部掛け合わせたもの(の nでの剰余)とするのです。 これにより、最大公約数を求めるのにかかる時間を 1/Mにすることができます。 ただし、全部掛け合わせると \bmod n 0になってしまう場合(≒複数の因数が M個の中で同時に発見された場合)の対処が追加で必要になります。

この手法は非常に大きな性能向上につながるため、高速化のためには事実上必須のテクニックです。 以下では、ポラードのロー法といったとき、暗黙でこのテクニックが適用されているものとします。

本題:私の提出コードの解説

ボトルネックの解析

ポラードのロー法のボトルネックは、疑似乱数を順次作っていく部分にあります。  x_iモンゴメリ表現から x_{i+1} = x_i^2 + C \bmod nモンゴメリ表現を求めるには、モンゴメリ乗算が1回必要です。 モンゴメリ乗算には、以下が必要です。

  • 64bit×64bit=128bitの乗算3回
    • フルの128bitが必要なのが1回
    • 下位64bitだけ得られれば良いのが1回
    • 上位64bitだけ得られれば良いのが1回

これらの乗算はいずれも、IntelSandy Bridge以降のプロセッサでは3cycleで実行可能です。 これらの乗算は直列に依存しているため、これだけで9cycleが必要です。 さらに一回の加減算が必須で、あわせて10cycleのレイテンシが避けがたいです。

差を求めて掛け合わせる部分でもモンゴメリ乗算が1回必要なので、合計で一つの乱数生成当たり6乗算必要(ブレントの改良方式の場合)ということになります。 乗算器はパイプライン化されているので、10cycleの間にスループットとして10回の乗算が可能であり、たった6回しか乗算を行わないのは非常にもったいないです。 このように、疑似乱数を順次作っていく部分がクリティカルパスとなっており、並列に実行できる命令が十分にないため、CPUはかなりの時間遊んでいます(とはいっても平均IPCが2を超えるような命令レベル並列性が豊富なプログラムなのですが、最近のCPUの最大IPCは4以上あるので、それには届かないということです)。

ボトルネックの解決策

乗算器が余っていることに着目し、並列に二つの疑似乱数系列を生成していくようにすることで、高速化が期待できます。 このようにすると、ループ中に12乗算含まれる形になります。 IntelのCPUは1cycleにスループットとして1回の乗算が可能であり、このループの逆数スループットは12cycleとなります。 これは疑似乱数生成のレイテンシは10cycleよりも長いので、理想的にはスループットに支配されたプログラムとなります。 ただし、逆数スループット12cycleはクリティカルパスレイテンシ10cycleとあまり差がないため、あまりにも雑な命令スケジューリングが行われるとレイテンシが見えるようになってしまう点に注意する必要があります。 普通のCPUの命令スケジューリングでは古い命令が優先的に選択されることを考えると、掛け合わせる部分の命令列よりも十分前に疑似乱数を生成する部分の命令列があれば良さそうです。 この目的のためには、ソフトウェアパイプライン化を施せばよいです。

実は、以下の条件を満たす場合はソフトウェアパイプライン化を施さなくてもそこそこの命令列が生成されていました。

  • clangでコンパイルする
  • 掛け合わせる部分のコードが実行されるのが、乱数生成直後ではない(例えば、差を計算するのに時間がかかるなど)

gccコンパイルすると遅いのはなんでだろうとずっと悩んでいたのですが、clangを使っていても差を計算する部分を軽量化(後述)したときにかえって遅くなったことから、命令スケジューリングに問題があることに気づくことができました。 ソフトウェアパイプライン化を施したコードは、上記条件を満たすときのclangの出力コードよりも5%程度高速であり、これが単独fastestを取る詰めの一手となったようです。

(あまり効果のない)軽量化

モンゴメリリダクションしない

普通の感覚では、差を掛け合わせた数と nの最大公約数を求めるとき、まずは差を掛け合わせた数のモンゴメリ表現qに対し元の表現Qを求め(この手順をモンゴメリリダクションといいます)、その後Qnの最大公約数を求めると思います。 しかし、qQ 2^{64}倍した数を \bmod nで考えているだけなので、nの奇数の約数aについて、Qaの倍数であればqaの倍数です。 よって、元の表現Qを求めなくても、qnの最大公約数を求めれば目的を達成できます。

正規化しない引き算

普通の感覚では、モンゴメリ表現xモンゴメリ表現zの差を求める際、x < z ? x - z + mod : x - zなどと求めるはずです。 しかし、モンゴメリ表現は [0, n)の範囲に正規化されている必要はないため、x < zかにかかわらず常にx - z + modとしても大丈夫です。 また、xはかなり長い間固定なので、x + modを計算しておけば、差を求めるのはsub命令一つでできるようになります。

なお、これにより掛け合わせる部分のコードが実行されるのが乱数生成直後になってしまいかえって遅くなってしまったことが、命令スケジューリングに問題があることに気づけたきっかけでした。

謎の数のモンゴメリ表現を使う

普通の感覚では、 x_{i+1} = x_i^2 + C \bmod nの計算にモンゴメリ乗算を使う場合、 Cモンゴメリ表現を求める必要があるように思われます。 しかし、 Cの具体的な値は重要ではなく、 Cには任意の整数が使えます。 よって、「モンゴメリ表現が1になるような謎の定数 Cを使った疑似乱数生成器」を使っているんだと考えることで、モンゴメリ表現に変換するコードを記述しなくて済むようになります。 この定数 C nに依存するためアルゴリズムが正しく動作するのか怪しいところがありますが、少なくとも 2^{40}までの合成数の非自明な因数を見つけることは可能なようです(Intel Core i9-12900Kにて24スレッド並列で1日くらいかけて確認)。

この技術は、あらゆるところに使われています。 乱数のシードZ1, Z2モンゴメリ表現である必要がありますが、いかにも普通の表現のような顔をして1などと書かれています。 掛け合わせていくときの初期値q1, q2は、普通は1のモンゴメリ表現を使いますが、ここでもわざと普通の表現と混同して1などと書いています。

これらは、プログラム中の生存変数を減らすことを目的に行っています。 ただ、頑張って減らしても結局gccが発生させるスピルを取り除くことはできませんでした(clangだと生存変数が多くてもスピルも無駄なmov命令もほぼなかったのに……)。 最初はgccの生成するコードが遅い原因をスピルなのではないかと疑っていた(命令スケジューリングに問題があることに気づいていなかった)ため、生存変数を減らす努力をしていたのですが、結構無意味な努力だったようです。

疑似乱数生成器を固定する

一般に一つの疑似乱数生成器だけでは素因数分解ができない合成数があるため、疑似乱数生成器を[n]( uint64_t x ) { return __uint128_t(x) * x % n + 1; }などとハードコーディングするのは正しくない(無限ループしうる)ことを以前に指摘しました(間違ったポラードのロー法の使い方 - よーる)。 しかし、今回は二つの疑似乱数生成器を並列で使っているため、このような問題は多分起こりません(確証なし)。 よって、プログラムの命令数を減らすために、疑似乱数生成器をハードコーディングすることができそうです。 プログラムの命令数を減らすことが高速化につながっているのかはよくわかっていません。

ポラードの方法(フロイドの循環検出法)でもブレントの方法でもない循環検出法

そんな大層なものではなく、ブレントの方法の亜種です。

周期を Pとして、循環を検出する( nとの最大公約数が1でない「差」を見てくる)のに必要な乱数生成回数は、フロイドの循環検出法では 3P回程度、ブレントの改良方式では 2P回程度です(正確には周期だけでなく循環に達するまでの反復回数(ρのしっぽの部分の長さ)を考えなければなりません。これに関するわかりやすい説明が前掲のポラード・ロー素因数分解法について #競技プログラミング - Qiitaに載っています)。 このことからブレントの改良方式が優れているとの考えが主流ですが、これは必ずしも正しくありません。 二つの手法における計算量を確認してみます。

フロイドの循環検出法では、一周当たり20cycle(乱数生成を直列に2回やるレイテンシ。もう一系列の乱数生成のレイテンシは並列実行で隠蔽可能)で乗算が12回(3乗算×乱数生成3回+3乗算×掛け合わせ1回)です。 循環検出までにかかる周回数は P回程度です。

ブレントの改良方式では、一周当たり10cycle(乱数生成を1回やるレイテンシ)で乗算が6回(3乗算×乱数生成1回+3乗算×掛け合わせ1回)です。 巡回検出までにかかる周回数は 2P回程度です。

このようにしてみると、どちらが特に優れているなどとは言えなさそうに見えてくると思います。

ブレントの改良はもう一つあり、差を取る乱数の片側(提出コードにおけるx側)を変更した直後は掛け合わせを省くという最適化です。 前掲のポラード・ロー素因数分解法について #競技プログラミング - Qiitaによれば、

もし r/2 < k < 3r/4 の間に約数が検出されるとすると、残りの 3r/4 ≦ k < r の間でも再度同じ約数が検出されるため

という説明がされています。 これはアルゴリズムの正当性(および反復回数を上から抑える議論)としては正しいものですが、近年のCPUにおいて高速化することの説明としては正しくなさそうです。 ブレントの改良方式では、一周当たり10cycleで乗算が6回なので、ここから掛け合わせを省いて乗算を3回に減らしても特に高速化しないからです(クリティカルパス上に無い命令が減ってもレイテンシに支配されているプログラムは高速化しません)。 逆にこの“スキップ”により、せっかく乱数を生成しているにもかかわらず因数を見つけるチャンスを失っていることになります。 この機会損失と最大公約数計算の手間のどちらが大きいかは簡単には分かりませんが、単純に「仕事量が減っているから高速化するはず」という雑な議論では済まなさそうです。

私の導入した二つの疑似乱数系列の並列生成を導入した場合、以下のようになります。

  • 掛け合わせを省かない場合、一周当たり12乗算なのでスループットで支配されて12cycleかかる。よって一周は12cycle
  • 掛け合わせを省く場合、一周当たり6乗算なのでスループットは十分で、レイテンシの10cycleがみえる。よって一周は10cycle

つまり、掛け合わせを省くことで多少は高速化するようです。 しかし、実験してみるとこのような先ほどよりも良い条件にも関わらず、素因数分解の速度の向上は見られませんでした。 やはり、“スキップ”中の因子発見チャンスを捨てている損失がある程度存在するようです。

フロイドの循環検出法を用いたものも実験してみましたが、ほぼ同じかやや遅いくらいの性能になりました。 原因は正確には調査できていないですが、どうも「全部掛け合わせると \bmod n 0になってしまう場合」という運の悪いパターンの発生率が上がっているようでした。 ソフトウェアパイプライン化がうまくいっていないなどの原因も考えられます。

以上の実験により、提出コードは“スキップ”を省いたブレントの改良方式を使っています。

パラメータの調整

差を512個掛け合わせてから最大公約数を求める

差を掛け合わせる個数 Mを大きくすることには、以下のような得失があります。

  • 最大公約数を求める回数が減るため、高速化する
  • 「全部掛け合わせると \bmod n 0になってしまう場合」の発生率が上がる(発生するとかなり低速になる)

実験を行った結果、 M=512程度を選ぶと実行時間が最小化されるという結論に達しました。

最悪実行時間に最適化

ブレントの改良方式では、差を取る乱数のうち固定側(提出コードにおけるx側)を更新するのは乱数生成の1, 2, 4, 8, ...回目にするのが普通です。 しかしこれは増やしていくという部分が重要なのであって、最初は1回で更新、などとする必要はありません。 逆にこのようにすると、頻繁に最大公約数を求める手順が挟まるため、大きな合成数素因数分解では不利そうです。 よって、提出コードでは、乱数生成の512, 1024, 2048, 4096, ...回目に更新することにしました。 これにより、小さい合成数に対する速度がかなり低下しますが、(少なくともこの問題に対しては)小さい合成数に対する速度が上がっても得にならないので、このような最適化を行いました。 そもそも、小さい合成数素因数分解を速く行いたいのなら、ロー法ではなく試し割り法の方が向いているという事実もあります。

その他

決定的素数判定法のコードに関する疑惑

本題と関係ありませんが、ミラーラビン素数判定法でif( n <= a ) { return true; }とやるのはまずいという話を見つけました(ミラーラビン素数判定法についてメモのコメント)。 これは一般には多分そうなのですが、提出コードで使われている64bit整数に対する決定的素数判定法(競技プログラミングでよく使われる方法)では問題ないはずです。 まず、n == 2の場合は素数、そうでなくn % 2 == 0の場合は素数ないと判定する部分は問題ありません。 その後、まずはa = 2でミラーラビン素数判定が行われます。 底2での強擬素数(ミラーラビン素数判定が誤判定する合成数)は最小でも2047なので、このテストが通っている時点でif( n <= 7 ) { return true; }if( n <= 61 ) { return true; }は正当です。 七つの底を使う方では、そもそもそちらが実行される条件がn >= 4'759'123'141なので、if( n <= 1'795'265'022 ) { return true; }などの条件が成立することは決してありません。

プログラムの実行時間の比率と乗算器利用効率

おおよそ以下のような感じでした。

  • 素数判定1.5%
  • 最大公約数を求める部分4.3%
  • 最大公約数を求め終わった時の分岐予測ミスペナルティ0.4%
  • 乱数生成部分93.1%
  • その他0.7%

乱数生成部分は一周あたり12cycleで512回繰り返されるので6144cycle、最大公約数の計算は200cycle程度なので、比率のオーダー感としてはあっています。 分岐予測ミスペナルティは20cycleくらいなので、最大公約数の計算にかかる時間との比率もあっています。 乱数生成部分で乗算器を使った回数をカウントし、それを実行にかかった時間で割って得た乗算器利用効率は、86.3%でした。 乱数生成部分が占める時間は93.1%なので、92.7%ほどは乗算器を有効利用できているということになりそうです。

なぜ二並列なのに二倍速くならないのか

まず、二並列にした時点でレイテンシ支配ではなくスループット支配に変わっています。 レイテンシは一周当たり10cycle、逆数スループットは12cycle(つまり一系列当たり6cycle)なので、この時点で1.67倍までしか速くならないことになります。

また、ロー法で循環を検出するまでにかかる平均周回数が半分になるわけではないのも原因のようです。 ロー法で循環を検出するまでにかかる周回数が指数分布に従っているとすれば、(仕事を並列化するのではなく)二系列同時にやって早く巡回を検出したほうを採用、という手順で平均周回数が半分になって、二倍高速化するはずです。 しかし、ロー法で巡回を検出するまでにかかる周回数は指数分布から外れているらしく、二系列同時にやっても平均周回数は半分にはならず、2/3程度にしかならないようです。 分布の標準偏差は平均値の半分程度(指数分布なら平均値と同じくらいのはず)であり、短い巡回が見つかる可能性がかなり低いことを示しています。 非巡回部分(ρのしっぽの部分)に原因があるのかもしれません。 また、歪度は1程度です(これは短い巡回には1という自明な下限があるにもかかわらず長い巡回には上限がないことが原因で、普通です。指数分布なら2あるはずなので、やはり短い巡回が少ないことを意味しているようです)。 尖度は1(四次中央モーメントで書くと4)程度であり、裾の重い分布であることを示しています。

これらから計算すると、提出コードは二並列にしているにもかかわらず実行時間は 2/3\times12/10 = 0.8倍程度、あるいは実行速度は1.25倍程度、にしかならないとわかります。

まとめ

  •  10^{18}までの整数(60bit整数)を素因数分解するプログラムを高速化し、Library Checkerにおいて最高速を達成した
  • ECM楕円曲線を用いた素因数分解法)を使うのではなく、ポラードのロー法をチューニングした、非常にわかりやすい実装
  • 乱数を連鎖的に生成する部分がレイテンシ支配になっていてもったいないことに着目
  • 二つの乱数系列を並列に生成することで高速化できた
  • レイテンシが見えてこないよう、ソフトウェアパイプラインを組む必要があった
  • ブレントによる改良型循環検出アルゴリズムは、別にそんなにすぐれているわけではない。また“スキップ”はやらない方がよさそうだった

*1: x_i^2の部分で正負の違いが無視されて同じ系列に合流するということは、巡回には正の数と負の数のどちらもが入ることはないということで、周期の最大は (a+1)/2な気がします。

*2:trailing zero countが3cycle、引き算命令とシフト命令が1cycleと仮定

Unityのgit管理をWSL2でやる

最近久しぶりにUnityでゲーム作りをしています。

共同で開発するためにgitを使うのですが、gitコマンドをWSL2から使用するといまいちな感じになったので、何とかする方法のメモです。

プロジェクトはWSL2の外に作る

UnityプロジェクトをWSL2の中(/home/lpha/以下など)に作るとうまくいかなかったので、Windows領域(/mnt/c/Users/lpha/Desktop/以下など)に置くことにしました。 いちいちcd /mnt/c/Users/lpha/Desktop/などと打ち込むのも面倒なので、cdsというエイリアスを作りました。

パーミッションの管理

WindowsにはUNIX系のようなファイルパーミッションの概念がないようで、この部分で問題が発生します。 新しく作ったファイルやUnityが自動生成したファイル(.metaファイルなど)のパーミッションは、VisualStudioやSourceTreeなどのWindows向けgitクライアントからコミットすると、644になります。 一方で、このリポジトリをWSL2のgitコマンドで見ると、コミット直後にもかかわらず変更があるように表示されます。 これは、Windows領域のファイルが755扱いで見えているからです。 つまり、gitの認識している状態(HEADではファイルのパーミッションは644)とWSL2の認識している状態(ファイルのパーミッションは755である)が異なるので、変更があるように表示されます。

また、git reset --hard HEADなどとしても、変更がある状態のままになってしまいます。 これは、WSL2はWindows領域のファイルのパーミッションを(デフォルトだと)記録せず、常に755であると認識するためです(記録されないので変更することもできません)。 これに対しては対処方法があります。 /etc/wsl.confに以下の記述を行うことで、Windows領域のファイルのパーミッションを記録することができるようになります(metadataの部分が重要)。

[automount]
options="metadata,umask=22,fmask=11"

なお、この記述を加えた後、WSL2を再起動する必要があります。 そのためには、コマンドプロンプトから以下のコマンドを実行します。

wsl --shutdown

パーミッションの管理2

このように設定しても、デフォルトのパーミッションが一貫していない点は変わりがありません。 なので、VisualStudioやSourceTreeなどからコミットした直後は、WSL2のgitコマンドで見ると変更がある状態になっています。 これは、Windows領域で新規作成されたファイルのパーミッションは何も記録されていないので、WSL2としては755であると判断するしかないためです。 しかし、git reset --hard HEADなどとすれば、ファイルのパーミッションが644に変更され、“正しい”状態になります。 これにより変更がある状態は解消され、以降のgitコマンドはいつも通り使える状態になります。

また、WSL2のgitコマンドでgit commitする場合は、ファイルのパーミッションが755になっていないか気を付けたほうが良いでしょう。 そうしないと、パーミッションが変更されるだけのコミットができたりして、コミット履歴が汚くなってしまいます。 これは多分pre-commitフックなどを使って対処するべき問題のような気がしますが、コミット作成はGUIツールの方が楽なのでgit commitコマンドを直に使用することがなく、放置しています。

メモ

どうもパーミッションを無視するという設定もできるようです。 ただ、私がWSL2からgitコマンドを使用したいという理由で共同開発者全員に強制するのも忍びなく、WSL2を使う側の負担のみで可能な方法を紹介しました。

なぜドルコスト平均法しないほうがいいか

新年になりました。

タイトルに関してですが、なんだか話があったので軽くまとめておきます。

要点としては、

  • 投資の期待値が標準偏差に比して大きい(インデクスファンドなどはこの条件を満たす)投資商品がある場合、
  • 余剰資金があるなら、
  • 何も考えずにその資金を即投入したほうが良い

というものです。

これに対して、「ドルコスト平均法は(高値で買ってしまう)リスクを削減しているから有効なはずだ」という反論がありますが、これが間違っていることを説明します。

なお、「投資の期待値が標準偏差に比して大きい」が満たされない場合、ドルコスト平均法は一定の合理性があります*1

直感的な理解

ドルコスト平均法を行うと、現金で持っている期間が存在するため、リターンを得る機会を逃します。 行いたかったことは「投資を行って、リターンを得る」だったはずなのに、「高値で買うリスクを回避」に目的が変わってしまっています。

その上、「高値で買うリスクを回避」と言っているにもかかわらず、値上がりしそうな投資商品を選んでいるため、待っている間に値上がりする蓋然性が高いです。 待てば待つほど、実のところ、高値で買うリスクを引き受けていることにすらなります。

統計的な理解

たしかにドルコスト平均法はリスクを削減します。 しかしその効果のほとんどは、「現金で持っている期間が長い」ことによります。 ドルコスト平均法では投資額面が平均で半分になるため、リスクを半分に下げる一方で、リターンも半分に下がってしまいます。

ドルコスト平均法特有の効果としては、

  • 投資商品の分散が大きくなるほど、高値で買うリスクを回避したことによる、追加のリスク削減が得られる
  • 一方で、投資商品の期待値が大きくなるほど、複利効果を逃したことによる、追加のリターン削減が発生する

のふたつがあげられます。

半額を投資する方法との比較

リスクを半分にしリターンを半分にするなら、余剰資金の半分を期初で一括で投資に回し、残りを現金で持っておく、という方法でも同じ効果が得られます。 これに対しドルコスト平均法は、

  • 分散がやや下がる
  • 期待リターンがやや下がる
  • あわせると、シャープレシオはわずかに上がる
  • リターンの中央値はほぼ同じ
  • リターンの第一四分位・第三四分位は下がる
  • 期初でのリスクが低い
  • 期末でのリスクが倍
    • 期末には全額を投資に回していることになるので

という違いがあります。

*1:正確にはこれは誤謬で、流動性が低くなっている時に取引をしている(安くなっている時には売り手が少なく、高くなっている時には買い手が少ない)ため売買が成立しにくい、という問題があります。

世界一簡単な物理レジスタ方式アウトオブオーダーCPUの作り方

この記事は自作CPU Advent Calendar 2023の24日目の記事です。

みなさんが使っているパソコンのCPUは、そのほとんどがアウトオブオーダー実行をしているはずです。 しかし、アウトオブオーダーCPUの作り方が書かれた教科書はあまりないような気がします。 トマスロのアルゴリズムというものが紹介されることもありますが、これは現代のアウトオブオーダー実行方式とはかなり違います。 この記事では、より洗練された現代的なアウトオブオーダー実行をするプロセッサを作ってみます。

全体的な構成

今回作るプロセッサは、FPGA上で動き、RV32IZmmul命令セットを動かすことができます。 東大のプロセッサ実験と同じルール(命令メモリ32KiB、データメモリ64KiB)で作ります(プロセッサ設計で一番難しいキャッシュを考えなくていいのですごく楽です)。

パイプラインは、以下のように設計しました。

  • 分岐予測ステージ
  • フェッチステージ
  • デコード+リネーム+ディスパッチステージ
  • スケジュールステージ
  • 発行ステージ
  • 実行ステージ
    • 普通の整数命令のレイテンシは1サイクル
    • 乗算命令のレイテンシは2サイクル
    • ロード命令のレイテンシは3サイクル
  • コミットステージ

今回は、作るのが大変なのでスーパースカラではないアウトオブオーダープロセッサを作ります。 この記事の作り方はスーパースカラであっても適用可能ですが、スーパースカラだといろいろ考えることが多くて大変です。

なお、スーパースカラではないのでIPCはどう頑張っても1.0を超えません。 普通にインオーダーのスーパースカラプロセッサを作ったほうが高いIPCを出せます。 性能が高いプロセッサを作ることが目的なら、こんな面倒なことをしてアウトオブオーダー実行を導入する必要はないです。 今回の話は、今後スーパースカラのアウトオブオーダープロセッサを作る土台作りということでやっています。

以下、各ステージを詳しく見ていきます。 なお、プロセッサが同時に取り扱える命令の数(リオーダーバッファサイズ)は16としました。物理レジスタ数は48です。

分岐予測ステージ

プロセッサの中で性能に最も影響を与えるステージです。

方向分岐予測器には、hashed perceptron分岐予測器を採用しました。 Hashed perceptron分岐予測器は、高性能分岐予測器の中では作るのが最も簡単です。 現在最強と考えられているのはTAGE分岐予測器ですが、hashed perceptron分岐予測器もそれとほとんど変わらないような性能を達成します。 TAGE分岐予測器は実装がとても大変なので、今回はhashed perceptron分岐予測器を使うことにします。

分岐先予測器は、分岐先バッファ(BTB)のみとしました。 置換とかのめんどくさいことを考えたくないので、32KiBの命令メモリに入る8192命令すべてに対して4Byteの飛び先を記録できるだけのBTB(32KiB)を用意しました。 こういう無駄遣いをすることで、アウトオブオーダープロセッサを簡単に作っていきます。 消費資源量の最適化は、正しく動いた後でやりましょう。

精度高く分岐予測するためには、方向分岐予測器で使うグローバル履歴をまだ結果がわかっていないにもかかわらず、予測に基づいて=投機的に更新する必要があります。 しかも、分岐予測ミスが起こった時にこの投機的更新を正しく取り消さないと、それ以降の分岐予測がめちゃくちゃになります。 これの対策は結構簡単で、グローバル履歴をリングバッファに持っておき、各命令でリングバッファのheadポインタを控えておけばよいです。 今回採用したhashed perceptron分岐予測器の最大分岐履歴長は90にしましたが、グローバル履歴のリングバッファサイズはきりよく128にしました。 これにより、38個までの投機的更新を行っても、それを取り消すことが可能です。 実際には高々16回しか投機的更新は行われないため、全然余裕です。

分岐予測ステージの結果は、

  • 次のサイクルの分岐予測の“種”として、分岐予測ステージへ
  • 次のサイクルのフェッチするアドレスとして、フェッチステージへ
  • 今フェッチした命令の「予測されたnext pc」(これは分岐を実行した時に分岐予測が正しいか判定するために必要)として、デコードステージへ

それぞれ送られます。

フェッチステージ

32KiBの命令メモリから4Byte命令をフェッチするだけで、特に工夫はありません。 命令キャッシュではないため、失敗はありません(常に命令をフェッチすることができます)。

デコード+リネーム+ディスパッチステージ

デコード

デコードするだけです。なにもむずかしいことはありません。

今回の作り方では、命令のレイテンシ情報、ソースレジスタを何個持つか、デスティネーションレジスタを持つか、だけをデコードしています。 発行ステージでデコードをする余裕があるため、ALUコードなどはそこでデコードをすることにして、4Byteの命令をそのまま後段に送っています。 このように作ることで、デバッグがやりやすくなっています。

リネーム

レジスタリネーミングを行います。つまり、命令に含まれる論理レジスタ番号を物理レジスタ番号に変換します。 アウトオブオーダー実行した結果を格納するために物理レジスタを使うのが、現代的なアウトオブオーダー実行方式です。

そもそも、アウトオブオーダー実行した結果をそのままアーキテクチャレジスタ(論理レジスタのこと)に書き込んでしまうと分岐予測ミスなどが起こった時に取り返しのつかないことになるので、アウトオブオーダー実行した結果はアーキテクチャレジスタとは別の場所に書いておく必要があります。 これの古典的な実装として、「アウトオブオーダー実行した結果はリオーダーバッファに書き込んでおき、分岐予測ミス等で取り消されないことが確定し次第順番にアーキテクチャレジスタに書き戻す」というものがあります。 しかし、この方式はデータの移動が発生してうれしくないです。 これより優れているのは、アーキテクチャレジスタとアウトオブオーダー実行の結果を書き込むバッファを区別しないで、一体の物理レジスタを使う、という方式です。 論理レジスタの値が物理レジスタのどこに入っているかは、そのポインタだけを記録しておけば、データ本体を動かす必要がなくなります。 「indirectは全てを解決する」というやつです。

ソース論理レジスタ番号を物理レジスタ番号に変換するのは簡単です。 その対応関係が記録された表(レジスタマップテーブル、RMT)を参照するだけです。

デスティネーション論理レジスタ番号は、未使用な物理レジスタ番号に変換される必要があります。 なぜなら、使用中の物理レジスタをデスティネーションとして割り当ててしまうと、他の命令が読もうとしていたデータが上書きされて別のデータに変わってしまい、プロセッサの動作がおかしくなる可能性があるからです。 ここで、「未使用」とか「使用中」とかは、「記録されているデータがもう二度と読まれない」「記録されているデータがまだ読まれる可能性がある」の意味です。

未使用な物理レジスタ番号は、フリーリストとよばれるリングバッファに記録しておくことにして、割り当てが必要になった場合はそこから一つ取り出します。 一般には未使用な物理レジスタが存在しない場合がありますが、今回は十分な物理レジスタを用意しました。 よって、この操作も失敗することがなく、簡単に実装できます。

後続の命令のリネームが正しくなるように、この命令に割り当てられた物理レジスタ番号をRMTに書き込んで更新します。 例えば、この命令のデスティネーション論理レジスタ番号が13だとして、それに物理レジスタ番号41が割り当てられたとします。 RMTにこのことを書いておけば、後続の命令でソース論理レジスタ番号が13の命令は物理レジスタ番号41を読みに行けばいい、と知ることができます。

この書き込みを行う前に、デスティネーション論理レジスタ番号に「今」対応する物理レジスタ番号も読み出しておきます。 これは、この命令がコミット(確定)された時に、永久に読み出されなくなる物理レジスタです。 よって、この命令をコミットするときに、その物理レジスタ番号をフリーリストに記録すれば、「未使用な物理レジスタ番号がフリーリストに入っている」状態を維持できます。 このようにして、物理レジスタガベージコレクション(といってもstd::unique_ptrと同様に自明ですが)が行われます。 物理レジスタ番号をフリーリストから一つ取得した命令は、そのうち物理レジスタ番号を必ず一つフリーリストへ返却するため、未使用な物理レジスタが枯渇することはありません。 RMTに書き込んでしまうと何を返却すべきかの情報が失われるので、この時点で読み出して保存しておく必要があります。

さて、アウトオブオーダープロセッサで最も面倒くさいのは、分岐予測ミスが起こった時にRMTをその分岐命令の時点の状況に戻さないといけない点です。 今回は贅沢に、すべての命令の時点でのRMTをまるまるコピーしておく(すべての命令でチェックポイントを取っておく)という方法を採用しました。 これはかなり贅沢で、RMTだけで〈6bit×32エントリ×16バージョン〉も必要になっています。 これは物理レジスタファイルに匹敵するサイズであり、しかも大量のコピーが発生しているので「実際のデータを移動させずポインタを変えるだけでいいから優れている」とか言ってたのを満たせていない気がしますが、気にしないことにします。

メモリ依存予測

本来であれば、リネーム(命令間のレジスタを介した依存の解析に相当)に加えて、メモリ依存予測(命令間のメモリを介した依存の解析)を行ったほうが良いです。 アウトオブオーダー実行をしている以上、プログラム順で後にあるロード命令が、プログラム順で前にあるストア命令よりも先に実行されることがあります。 この時、アドレスが一致していれば実行結果がおかしくなっている可能性が大です(メモリ順序違反)。 これを検出する機構は当然搭載されるため、間違った実行結果にはなりません(そういうことが起こったら分岐予測ミス同様にフラッシュしてやりなおします)。 とはいえ、頻繁にフラッシュが発生すると性能に影響があるので、「このロード命令はこのストア命令の後にやった方がいい」みたいな予測(これをメモリ依存予測と言います)をするわけです。

ただ、メモリ依存予測を実装するのはすごく大変なので、今回は実装しないことにしました。 実際の所、まともに最適化が行われていれば、ストアしてすぐロードすることはほぼありません。 したがって、命令ウィンドウ(命令を並べ替える範囲)が小さければ、メモリ順序違反はほぼ発生しないため、メモリ依存予測を行わなくてもそれほど性能が落ちることはありません。

また、メモリ順序違反を検出する機構であるロードストアキューのエントリの確保もこの辺のステージで行うのが普通ですが、今回は十分な数(16)のエントリを用意しました。 よって、確保が失敗することはなく、簡単です。

ディスパッチ

スケジューラ(どの命令が実行可能であるかを判断し、その中から適切なものを選ぶ機構)に命令を登録します。 また、スケジューリングに関係ない情報(命令のPCや即値など)を脇によけておきます(payload RAMに格納しておきます)。 これは、スケジュールステージは一番複雑なステージであり、なるべく小さく作りたいからです。

普通のスケジューラの作り方をすると、「この命令が実行された(のでこの命令に依存している命令は実行可能になった)」という信号は1サイクルしか流れません。 なので、

  • ディスパッチをする前に「すでにこの命令は実行されている(のでこの命令に依存している命令はとっくに実行可能である)」という情報を手に入れる必要がある
    • しかもこの情報を適切に更新する必要がある
  • ディスパッチをするその1サイクルに流れる「この命令が実行された(のでこの命令に依存している命令は実行可能になった)」という情報を逃さないようにする必要がある

という二点をちゃんと作らないといけないので大変です。

今回は、「どの物理レジスタが読み出せる状態になっているかを管理する」という贅沢な方式でスケジューラを作ることで、ディスパッチを大幅に簡単に作れるようにしました。 また、スケジューラに入れられる命令数を16と十分に大きくとったため、ディスパッチが失敗することもなく、やはり簡単に作ることができます。

このステージの結果(高々二つのソース物理レジスタ番号と高々一つのデスティネーション物理レジスタ番号と命令の実行レイテンシ情報)はスケジューラに書き込みます。 また、デスティネーション物理レジスタについて、「この物理レジスタはまだ読み出せる状態ではない」と記録します(表に0を書き込みます)。

スケジュールステージ

スケジューラ内にある高々16命令のうち、どの命令が実行可能であるか判断し、実行可能な命令の中から最も古い命令を一つ選びます。

スケジューラ内にある16の命令は各々、エントリのvalidビットに加え「自分のソース物理レジスタは読み出せる状態か」を二つ確認し、それらのandを取ることで実行可能状態(ready)であるかを判定します。 スケジューラはリングバッファになっていて、古い命令から順番に並んでいます。 よって16本のready信号を確認して、その中で最も古いものがどれかを調べることはそれほど難しくありません。 これにより、実行可能である最も古い命令を一つ選ぶことができます。

選ばれた命令は、自身のエントリのvalidビットを0にすることで、再び選択されることを防ぎます。 その後、選ばれた命令のレイテンシー1サイクルだけ待ってから(つまり、レイテンシが1サイクルの命令なら選ばれたサイクルの終了時に)、選ばれた命令のデスティネーション物理レジスタについて、「この物理レジスタは読み出せる状態である」と記録します(表に1を書き込みます)。

スケジューラが管理してる「物理レジスタは読み出せる状態である」は実際には2サイクル(1サイクル+発行レイテンシ分)未来の状態を表しています。 普通はこれはレイテンシ予測に基づくもので、未来の状態が本当にそうなるかはわかりません。 未来予測が正しくなくなる典型的なパターンはキャッシュミスした場合です。 キャッシュミスすればデータは得られないので、「物理レジスタは読み出せる状態である」というのは間違いになります(レイテンシ予測ミス)。 間違った未来予測に基づいてスケジュールすると、正しいデータを使って実行できないので、プロセッサの動作がおかしくなります。 したがって、未来予測が間違った場合の対処を考える必要があります。 これは、アウトオブオーダープロセッサを作るうえで、難しい問題です(リザベーションステーションを使った方式はこの問題が起こらないという利点が一応あります。また、レイテンシが不確かなものについては実際にデータが得られるまで表を更新しない(投機的発行しない)という方式もあります)。

しかし、今回はデータキャッシュがなくてレイテンシ予測は完全に当たるので、簡単に作ることができます。

「物理レジスタは読み出せる状態である」の表は48個ある物理レジスタに対して1bitの情報を管理しているだけですが、ポート数がすごく多いので贅沢なつくりになっています(とはいえFPGAを前提にするとこれより小さく作る方法はそう簡単にはなさそうです。この作り方は結構いい方法なのでは?)。 この表は、スケジューラ内の16個の命令が2つのソースについて確認するので、32readになっています。 しかも、ディスパッチの1writeに加えて、レイテンシが1cycle/2cycle/3cycleの三パターンがあるので3writeが必要で、合計で4writeになります。 つまり、この表だけで〈1bit×48エントリ×32read×4write〉が必要になっています。

ここで選ばれた命令(実際にはその命令の情報が格納されているpayload RAMのインデクス番号)は、発行ステージへ送られます。

発行ステージ

Payload RAMを読み出し、PCと命令語などを得ます。 命令語はここでデコードします。 また、物理レジスタからの読み出しを行います(分散RAMで作って非同期読み出しします)。 さらに、必要であればフォワーディングを受け取ります。 そして、それらのうち適切な方をパイプラインレジスタに書き込みます。

実行ステージ

整数命令

命令を実行するだけです。 命令の実行結果は、即フォワーディングネットワークに流します。

一方、レジスタファイルへの書き込みはまだ行いません。 これは、レジスタファイルの書き込みポートを節約するためです(後から考えると、これは「資源を贅沢に使って簡単に作る」に反していました……)。 整数命令の実行レイテンシは1サイクル、最長の実行レイテンシを持つのはロード命令の3サイクル、ということで、差分の2サイクルは待ちます。 これは普通の五段パイプラインと同じ設計です。 レジスタファイルに書き込まれるまでの間は、フォワーディング経由で値を受け取れるので、これが問題になることはありません。

また、レジスタファイルに書き込むのと同時に、リオーダーバッファに実行完了フラグを書き込みます。 リオーダーバッファといっても、実行結果データを書き込むタイプではありません。

もし分岐命令であれば、実行の結果得られたnext pcが予測結果と一致しているかを確認します。 予測結果と一致していてもしていなくても、分岐予測器を学習します。 また、予測結果と一致していなければ、分岐予測ミスです。 例外の発生するサイクルをロードストア命令と合わせるため、第二サイクルで分岐予測ミス信号を送出します。 これは、複数の実行ユニットから出てきた例外の中で最古のものを選ぶ、という面倒事を回避するためです。 さすがに次のサイクルから正しい命令をフェッチするのはレイテンシ的に厳しいため、このサイクルではmispredフリップフロップに1を書き込むのにとどめ、次のサイクルから分岐予測ミスの対処を開始します。

乗算命令

2サイクルかけて乗算を行うだけです。 乗算結果は即フォワーディングネットワークに流しますが、レジスタファイル・リオーダーバッファへの書き込みは、1サイクル待ちます。 ロード命令に合わせて3サイクルにしておけばこの辺やスケジューラを多少簡単にできたので、レイテンシを2サイクルにしたのは失敗でした(乗算の結果をまた乗算に使うということはあまり発生しないので、レイテンシが長くてもIPCへの影響は軽微です)。

ロード命令

第一サイクルではメモリアドレスを計算します。TLBがないので単純です。

第二サイクルではメモリアドレスを使って、データメモリから読み出します。キャッシュがないので単純です。 さらに、メモリアドレスを使って、ロードストアキューを連想検索します。 もし、ロードストアキューに当該メモリアドレスに対するストアが記録されていれば、その中で最新のストア命令の書いた値を読み出します(ストアフォワーディング)。 また、ロードストアキューに、アドレスを記録します。

第三サイクルでは、ロードストアキューを連想検索して値が得られたならそれを、得られなかったならデータメモリから読み出した値を、零拡張または符号拡張します。 その結果は、フォワーディングネットワークに流します。 また、第三サイクルが終わる瞬間にレジスタファイル・リオーダーバッファに書き込みます。 これにより、発行ステージでのレジスタ値の選択が「レジスタファイルから読み出した値」「今完了した整数命令の結果」「今完了した乗算命令の結果/1サイクル前に完了した整数命令の結果」「今完了したロード命令の結果/1サイクル前に完了した乗算命令の結果/2サイクル前に完了した整数命令の結果」の4択になります。

※ストアフォワーディングはパーシャルメモリアクセス(SB命令で1byteストアしたのにLD命令で4byteロードするなど)を考えると面倒です。メモリアドレスの下位2bitごとにロードストアキューを作ることで対処しています。

ストア命令

第一サイクルではメモリアドレスを計算します。TLBがないので単純です。

第二サイクルでは、メモリアドレスとストアデータをロードストアキューに書き込みます。 さらに、メモリアドレスを使って、ロードストアキューを連想検索します。 もし、ロードストアキューに当該メモリアドレスに対するロードが記録されていれば、メモリ順序違反が起こったことを示しています。 その中で最古のロード命令を探し出し、その命令からフェッチしなおします(リフェッチ)。 このサイクルではmispredフリップフロップに1を書き込むにとどめ、次のサイクルからメモリ順序違反の対処を始めます。

第三サイクルは何もする必要がありません。

コミットステージ

アーキテクチャステートを不可逆的に更新するのがコミットステージです。 リオーダーバッファを確認し、最古の命令に対応するエントリに実行完了フラグが立っていればコミットを行います。 ストア命令であれば、ストアキューに書き込まれている値を読み出し、データメモリに書き込みます。 キャッシュがないので、これは確実に成功します。 ストア以外の命令であれば、何もする必要はありません。

その後、リオーダーバッファのheadポインタ*1 を1インクリメントします。 このheadポインタは、以下の情報を得るために使われます。

  • リオーダーバッファ中の最古の命令はどれか
  • スケジューラ内に入っている命令の中で「新しい命令」と「古い命令」の境界はどこか

また、分岐予測するごとに1インクリメントされるtailポインタと組み合わせて、以下の判定に使われます。

  • リオーダーバッファに空きがあるか
    • 分岐予測を進めていいかの判定で必要
  • ロードストアキューの中で有効な範囲はどこか
    • メモリ順序違反が起こったかの判定で必要

本当はここで、記録されている物理レジスタ番号をフリーリストに返却する必要があるはずです。 しかし、実は逆の発想で「リオーダーバッファの(headではなく)tailを進めるときに、物理レジスタ番号をフリーリストに返却」としても大丈夫です。 こうすると物理レジスタ番号がフリーリストに返却されるのが遅くなってうれしくないですが(というか物理レジスタ不足でストールが発生してtailが進められなくなればデッドロックしてしまう)、今回は十分な物理レジスタがあるため問題ありません。 パイプラインフラッシュが起こった時のことを考えるとこうしたほうがやりやすいため、この方式を採用しました。

パイプラインフラッシュの実装

さて、ここまででパイプラインの“正常”動作を書いてきましたが、実際には分岐予測ミスやメモリ順序違反などの投機ミスが起きるので、パイプラインフラッシュが発生します。 この時にプロセッサの状態を正しく戻さないと正しく動かないので、頑張って作ります。 といっても、今まで贅沢に作ってきたのはパイプラインフラッシュが起こった時にほとんど何もしなくていいようにするためなので、実装はそんなに難しくありません。

分岐予測ステージ

  • グローバル履歴のリングバッファのポインタを戻す
    • 一応、これをちゃんとやらなくてもプロセッサは正しく動くけれど、分岐予測がめちゃくちゃになる
  • 次のサイクルの分岐予測の“種”となるPCをセットしなおす

フェッチステージ

  • 正しいPCを受け取る(まだフェッチはしない)

デコード・リネーム・ディスパッチステージ

  • RMTの状態をチェックポイントから復元する(ポインタを正しい値に修正するだけ)
  • 今デコード中の命令は捨てる

スケジュールステージ

  • フラッシュの原因となった命令以降の命令のエントリについて、有効な命令が入っているかを修正する(validビットを0にする)
  • 物理レジスタのreadyの状態は、戻さなくてよい(ディスパッチ時に正しい値に修正されるので)

発行ステージ

  • 何もしなくてよい

実行ステージ

  • リオーダーバッファのheadとtailの範囲に無い命令であれば
    • mispredフリップフロップに書き込まないようにする
      • これは必須
    • レジスタファイルへの書き込みを行わないようにする
    • フォワーディングネットワークに値を流すのもやめる
      • この二つは今回のパイプライン設計だとやらなくても大丈夫な気がする

コミットステージ

  • 何もしなくてよい(もともとリオーダーバッファのtailより先をコミットしないので、tailを戻しておけば大丈夫)

物理レジスタガベージコレクション

分岐予測が行われて命令がリオーダーバッファに格納される(本当はまだ命令がフェッチされていないのでエントリが確保されるだけ)たびに、以下のどちらかを行えばよいです。

  • 「フラッシュされた」ビットが立っていれば、その命令が上書きするリオーダーバッファエントリに記録されている命令のデスティネーション物理レジスタ番号をフリーリストに返却(投機ミスからの回復動作)
  • 「フラッシュされた」ビットが立っていなければ、その命令が上書きするリオーダーバッファエントリに記録されている命令のデスティネーション論理レジスタにもともと対応していた物理レジスタ番号をフリーリストに返却(通常動作)

こうすることで、フリーリストの書き込みポート数を1のままにすることができます。 というかよく考えると、リオーダーバッファに記録されている情報(フラッシュビット+旧物理レジスタ番号+新物理レジスタ番号)がフリーリストの情報を特定するのに十分なものになっているので、フリーリストをわざわざ実装する必要はなかったようです……。

アウトオブオーダープロセッサの実行の様子

スーパースカラでないので、レイテンシが1サイクルの命令が並んでいるところはインオーダープロセッサと同じ感じになります。 CoreMarkというベンチマークを動かしてみたところ、以下の二か所ではアウトオブオーダー実行がうまく働いていそうでした。 パイプラインは、Konataというパイプライン可視化ツールで可視化しました。

図1: 乗算を行っている部分のパイプラインチャート。暇なサイクルでmove命令が実行できている(分岐予測ミスが起こるかもしれないのに先に実行できる!)。でもこの命令よく見ると無意味では……?
図2: ロードの連鎖がある部分のパイプラインチャート。5命令に6サイクルかかる(これは最古の命令を優先するという戦略がいまいちなため)のでリオーダーバッファ不足で性能が下がっているわけではない

まとめ

キャッシュなどの複雑化の要因をなくし、それでも複雑な部分は大量の資源を使って解決することで、簡単にアウトオブオーダープロセッサを作ることができることがわかりました。 来年はちゃんと省資源化(とできればスーパースカラ化)してみたいと思います。

*1:いろいろな流儀があるようですが、一番古い命令が「待ち行列の先頭」なのでhead、一番新しい命令が「待ち行列の末尾」なのでtail、という流儀を採用しています。パイプラインチャート的にも上がheadになっていい感じです。でも、headの方がなんか最新感があるし、動物は普通は頭がある方向に進む(蛇とか魚とか)のでheadがキューの先頭(のびる方向)であり一番新しい命令を指すべき、という意見も否定しがたいです。

0.1+0.2≠0.3について

浮動小数点数を扱う時に注意しなければいけない例として、0.1 + 0.20.3にならない、というものがあります。 具体的には、0.1 + 0.20.30000000000000004になってしまい、0.3にならないプログラミング言語が多数存在します。 この問題について、以下の記事では十進数から二進数に変換するときの誤差が原因であるとしています(この説は他にも多くのウェブサイトで唱えられています)。

qiita.com

たしかに十進数から二進数に変換するときの誤差は一因なのですが、より大きな寄与を持っているのは足し算の丸め誤差です。 そのことを確認していきたいと思います。

(わかっている読者への注:以下では、浮動小数点数と言ったらIEEE754で規定されたbinary64(二進倍精度浮動小数点数)を指すこととします)

浮動小数点数の基礎

基本的には、分子が 2^{52}以上 2^{53}未満の整数で、分母が二の整数乗であれば、浮動小数点数として表せます。 この形式で表せない数は、浮動小数点数として表せる数で近似します。 この近似を行うことを、丸めと言います。 普通は、最も近い浮動小数点数で近似します(最近接丸め。端数が0.5ぴったりの時のタイブレークの方式の違いに四捨五入や偶数丸めがあります)。 特別な用途(例えばお金が虚空から生まれてはいけないとか*1)では、「真の値を超えない範囲で」最も近い浮動小数点数に丸める(下向き丸め)などを使うこともあります。逆に「真の値を下回らない範囲で」とする方法もあります(上向き丸め)。

ULPはunit of last placeの略で、上記のような分数として表したときの、分子の1の違いを単位としたものです。 誤差は0.38ULP、のような形で誤差の解析によく使います。

計算の過程

0.1を浮動小数点数に変換する

0.1に最も近い浮動小数点数は、 \frac{7205759403792794}{72057594037927936} = \frac{\lceil0.1\times2^{56}\rceil}{2^{56}}です(十六進浮動小数点数リテラル表記で0x1.999999999999ap-4)。 この時点で、分子が7205759403792793.6であるべきところが7205759403792794になってしまっているので、+0.4ULPの誤差を生じています。

0.2を浮動小数点数に変換する

0.2に最も近い浮動小数点数は、 \frac{7205759403792794}{36028797018963968} = \frac{\lceil0.2\times2^{55}\rceil}{2^{55}}です(十六進浮動小数点数リテラル表記で0x1.999999999999ap-3)。 この時点で、分子が7205759403792793.6であるべきところが7205759403792794になってしまっているので、+0.4ULPの誤差を生じています。

0.3を浮動小数点数に変換する

0.3に最も近い浮動小数点数は、 \frac{5404319552844595}{18014398509481984} = \frac{\lfloor0.3\times2^{54}\rfloor}{2^{54}}です(十六進浮動小数点数リテラル表記で0x1.3333333333333p-2)。 この時点で、分子が5404319552844595.2であるべきところが5404319552844595になってしまっているので、-0.2ULPの誤差を生じています。

0.1と0.2を足す

丸める前の加算の結果は、 \frac{5404319552844595.5}{36028797018963968}です。 ここで、分子は 2^{52}以上 2^{53}未満ですが、整数になっていないので丸めを行う必要があります。 四捨五入なり偶数丸めなりを適用すると、 \frac{5404319552844596}{36028797018963968} = \frac{\lceil0.3\times2^{54}\rceil}{2^{54}}を得ます。

これは、0.3に最も近い浮動小数点数である \frac{5404319552844595}{18014398509481984} = \frac{\lfloor0.3\times2^{54}\rfloor}{2^{54}}と分子が1だけ異なります。

誤差の解析

0.1の指数部は0.3の指数部に対して2だけ小さいので、変換の誤差+0.4ULPが最終結果に与える影響は1/4となって+0.1ULPです。 0.2の指数部は0.3の指数部に対して1だけ小さいので、変換の誤差+0.4ULPが最終結果に与える影響は1/2となって+0.2ULPです。 加算の丸め誤差が最終結果に与える影響は+0.5ULPです。

よって、合計で+0.8ULPの誤差が生じます。

このように、この誤差の主要因は加算の丸め誤差です。

いつでもfaithfulになるか

Faithful rounding(信頼可能丸め*2)とは、真の値を上向き丸めして得られる値か真の値を下向き丸めして得られる値の、どちらかが返ってくるような丸めのことです。 言い換えると、真の値を挟む二つの浮動小数点数のうちのどちらかが返ってくるということです。 「最近接の浮動小数点数に丸めてね」というのは要求が厳しすぎるので、少し要件を緩和したということです。

先ほどの例では、真の値0.3に対し、最も近い浮動小数点数 \frac{\lfloor0.3\times2^{54}\rfloor}{2^{54}}を計算結果として得ることができませんでした。 しかし、得られた結果 \frac{\lceil0.3\times2^{54}\rceil}{2^{54}}は真の値を上向き丸めして得られる値です。 したがって、先ほどの計算は、最近接丸めではありませんでしたが信頼可能丸めではあったようです。

さて、十進小数同士の足し算結果が必ずしも最近接丸めではないということがわかりましたが、要件を信頼可能丸めに緩めれば満たされるでしょうか?

答えは否です。 まず、桁落ちが生じる場合は明らかにかなりの誤差が生じ得ます。 桁落ちが生じない場合に限っても、信頼可能丸めですらなくなる例は存在します。

0.17 + 0.28がその一例です。

0.17を浮動小数点数に変換すると、 \frac{6124895493223875}{36028797018963968} = \frac{\lceil0.17\times2^{55}\rceil}{2^{55}}となり、+0.44ULPの誤差が生じます。 0.28を浮動小数点数に変換すると、 \frac{5044031582654956}{18014398509481984} = \frac{\lceil0.28\times2^{54}\rceil}{2^{54}}となり、+0.48ULPの誤差が生じます。 真の値である0.45を浮動小数点数に変換すると、 \frac{8106479329266893}{18014398509481984} = \frac{\lceil0.45\times2^{54}\rceil}{2^{54}}となり、+0.2ULPの誤差が生じます。

二つの浮動小数点数が表す値を加算すると、 \frac{8106479329266893.5}{18014398509481984}となり、四捨五入なり偶数丸めなりをすると \frac{8106479329266894}{18014398509481984}となります。 これは真の値 \frac{8106479329266892.8}{18014398509481984} =  \frac{0.45\times2^{54}}{2^{54}}を挟む二つの浮動小数点数ですらありません。 つまり、加算を実行するだけで、信頼可能丸めにすらならないケースがあることがわかりました。

このケースの誤差は+1.2ULPです。 内訳としては、0.17を変換した時の誤差+0.44ULPが最終結果に与える影響は1/2されて+0.22ULP、0.28を変換した時の誤差+0.48ULPが最終結果に与える影響は+0.48ULP、加算の丸め誤差が0.5ULP、で合計+1.2ULPです。


このような場合、最大の誤差は1.25ULPとなります。 まず、変換の誤差は最大で0.5ULPです。 加算のオペランドの変換誤差が最終結果に与える影響は、最大で0.5ULPですが、両オペランドともそうなることはありません。 なぜなら、そのような場合は繰上りが発生して、変換誤差が最終結果に与える影響が1/2になってしまうからです。 よって最悪ケースでは、片方のオペランドの変換誤差が最終誤差に与える影響が0.5ULP、もう片方のオペランドの変換誤差が最終結果に与える影響が0.25ULPとなります。 これに加えて加算の丸め誤差が最大で0.5ULPなので、合計で最終結果に乗る誤差は最大で1.25ULPであるということになります。

最大誤差が1.0ULP未満であれば信頼可能丸めになりますが、1.0ULPを超えているため、信頼可能丸めにはなりません。 ただし、最大誤差が1.5ULP未満なので、「真の値を最近接丸めした値か、またはその前後の浮動小数点数」が計算結果となることは保証されます*3

*1:言うまでもなく、お金を扱う時に浮動小数点数を使うべきではないですが、例として

*2:忠実丸め、隣接浮動小数点数丸め、F丸め、などの訳も存在

*3:指数部の境界付近ではULPを使った議論は怪しいです。反例があったら教えてください。

gccが呼び出し規約に従わないコードを出力する例

アセンブリプログラムを書いていたら、遭遇したのでメモです。

さすがにプログラム全てをアセンブリ言語で書くのは大変なので、性能上重要ではないところはコンパイラのコード生成に任せて、性能上重要なところだけアセンブリプログラムで書くということはあると思います。 また、一からアセンブリ言語で書いたり全部イントリンシックで書くのも大変なので、コンパイラが出したコードをベースに改造することもあると思います。 この時、コンパイラが出すコードは呼び出し規約に沿っているので、性能上重要な関数を呼び出し規約に沿う範囲で自由に書き換えられると素朴に思っていたのですが、どうもそうではないということがわかりました。

問題を引き起こすコード

__attribute__((noinline)) char* id( char* x ) { return x; }

int main( int argc, char* argv[] ) {
    for( int i = 0; i < argc; ++i ) {
        argv[i] = id(argv[i]);
    }
}

gcc -O2コンパイルした場合、以下のコードになります(Compiler Explorerで確認したところ、gcc5~gcc13でどれもほぼ同じコードになりました)。

id:
        mov     rax, rdi
        ret
main:
        test    edi, edi
        jle     .L4
        movsx   rdi, edi
        lea     rdx, [rsi+rdi*8]
.L5:
        mov     rdi, QWORD PTR [rsi]
        add     rsi, 8
        call    id
        mov     QWORD PTR [rsi-8], rax
        cmp     rsi, rdx
        jne     .L5
.L4:
        xor     eax, eax
        ret

どこが問題か

このコードは、id関数を呼び出してもrsirdxが変わらないことを前提としています。 rsirdxは呼び出し側保存レジスタであり、id関数の中で書き換えても良いレジスタです。 つまり、このコードは呼び出し規約に沿っていないことになります。 おそらく、gccid関数でrsirdxが変わらないことを見抜くことで、このような積極的な最適化*1を行っているのだと思います。

こういうことをされたので、アセンブリファイルを勝手に書き換えたときにおかしくなったという話でした。

*1:ちゃんと呼び出し規約に従おうとすると、呼び出し先保存レジスタに値を入れる必要があって、その呼び出し先保存レジスタの確保のためにスピルが発生します。

x86のアドレッシングモードと使えないレジスタについて

x86CISC命令セットであり、RISC命令セットと異なりメモリ上のデータを命令のオペランドに取ることができます。 x86レジスタは(特に32ビットと64ビットは)おおよそ汎用ですが、一部のレジスタはアドレッシングに使えない制限があります。 よくわからなくなるのでまとめてみました。

16ビット

レジスタは全く汎用ではありません。基本的にはAXCXDXを算術演算に使い、BXSIDIをアドレス計算に使う、ということが想定されているようです。 基本的には、以下のアドレッシングモードが使えます。

  • [addr + offset]
    • addrは、以下の八種類
      • BX
      • BP
      • SI
      • DI
      • BX+SI
      • BX+DI
      • BP+SI
      • BP+DI
    • offsetは、以下の三種類
      • (何も足さない)
      • imm8
      • imm16
  • [imm16]

ただし、このうち[BP](つまり、addrBPで、offsetが「何も足さない」)は、そのエンコーディングが絶対16ビットアドレス指定に奪われているので、使えません。 もちろん、[BP+0]とすることで、1バイト増えてしまいますがエンコーディングすることは可能です。 なお、BPをフレームポインタとして使っている場合、[BP]に入っているのはリターンアドレスのはずなので、普通は使わないと思います。

32ビット

レジスタはほぼ汎用になりました。 基本的には、以下のアドレッシングモードが使えます。

  • [addr + offset]
    • addrは、以下の七種類
      • EAX
      • EBX
      • ECX
      • EDX
      • ESI
      • EDI
      • EBP
    • offsetは、以下の三種類
      • (何も足さない)
      • imm8
      • imm32
  • [imm32]
  • [addr + index * scale + offset](SIBバイト使用)
    • addrは、以下の八種類
      • EAX
      • EBX
      • ECX
      • EDX
      • ESI
      • EDI
      • ESP
      • EBP
    • indexは、以下の八種類
      • EAX
      • EBX
      • ECX
      • EDX
      • ESI
      • EDI
      • (零)
      • EBP
    • scaleは、以下の四種類
      • 1
      • 2
      • 4
      • 8
    • offsetは、以下の三種類
      • (何も足さない)
      • imm8
      • imm32
  • [imm32 + index * scale](SIBバイト使用)

16ビットの時からの伝統で、[EBP](つまり、addrEBPで、offsetが「何も足さない」)は、そのエンコーディングが絶対32ビットアドレス指定に奪われているので、使えません。 16ビットの時同様、[EBP+0]とすることで、1バイト増えてしまいますがエンコーディングすることが可能です。

SIBバイト利用アドレッシングでは、[EBP + index * scale](つまり、addrEBPで、offsetが「何も足さない」)は、そのエンコーディングが絶対32ビットベースアドレス指定に奪われているので、使えません。 同様に、[EBP + index * scale + 0]とすることで、1バイト増えてしまいますがエンコーディングすることが可能です。 indexEBPではなくscale1の場合は、逆順にすればバイト数増加はありません。

八つの汎用レジスタの一員であるにもかかわらず、[ESP + offset]は三種類とも、そのエンコーディングがSIBバイト利用アドレッシングに奪われているため、使えません。 これをやりたい場合は、SIBバイト利用アドレッシングの中にある[ESP + (零)* scale + offset]を使うことで、1バイト増えてしまいますがエンコーディングできます。 scaleは、1を設定するのが一般的です。

ESPは、そのエンコーディングが零レジスタ指定に奪われているので、indexとして指定することができません。 とはいえ、ESPを2倍した値というのが意味を持つ局面はないはずなので、通常は問題となりません。 なお、零レジスタを2倍・4倍・8倍する、というのも意味がないものであり、エンコーディングに無駄があります。

64ビット

汎用レジスタは16本に増えました。 アドレッシングモードは32ビットとほぼ変わりませんが、以下の三点に注意する必要があります。

  • [imm32]がなくなり、代わりに[RIP + imm32]ができるようになった
    • [imm32 + index * scale]は残存
  • RSPの“裏レジスタ”であるR12は、ModR/MでRSPが使えないことを引き継いでいる
    • indexRSPが使えない制限は引き継いでいない
  • RBPの“裏レジスタ”であるR13は、(ModR/Mとbaseで)RBPでできないことを引き継いでいる
    • indexではRBPは自由に使える

つまり、直観的でない点は、

  • [RIP + index * scale + offset]はできない
  • [R12][R12 + (零)*1]で表現されるので1バイト長い
  • [R12 + offset][R12 + (零)*1 + offset]で表現されるので1バイト長い
  • [R13][R13 + 0]で表現されるので1バイト長い
  • [R13 + index * scale][R13 + index * scale + 0]で表現されるので1バイト長い
    • scale1の場合は、逆順にすると短くできることがある

また、エンコーディングの無駄な点(rex.b01かで意味が変化しない)として、

が挙げられます(まぁそもそもrex.xはSIBバイトを使わない場合、rex.bはメモリオペランドを取らない場合、rex.wはバイト命令の場合、は無意味なのですが)。

まとめ

  • 16ビットの時は、アドレッシングに使えないレジスタが存在した
  • 32ビット以降では、基本的にすべてのレジスタがアドレッシングに使える
  • しかし、特定のレジスタを使うアドレッシングは、命令長が伸びることがある
    • [BP][EBP][RBP]は1バイトのびる(しかし、フレームポインタとして使っているならばそこに入っているのはリターンアドレスなので、あまり使われない)
    • [EBP + index * scale][RBP + index * scale]は1バイトのびる(しかし、フレームポインタとして使っていて、indexが素直な配列インデクスであれば、このパターンは出てこないはず)
    • [ESP][ESP + offset][RSP][RSP + offset]のようにESPRSPがmodR/Mで指定されるとベースだと1バイトのびる(これはよく出てくる)
    • ESPRSPindexには使えない(使われることはないはず)
  • 64ビットでは、“裏レジスタ”の関係にあるレジスタEBPESPのModR/Mとbaseにおける特性を引き継いでいる
    • indexでは制約はない
    • [R12][R12 + offset][R13][R13 + index * scale]は1バイトのびる