高速剰余算 div2by1 実装してみた
div2by1というアルゴリズムがある -> https://gmplib.org/~tege/division-paper.pdf
これはBarrett reductionやMontgomery乗算と違い、
- (k, k) -> 2k-bitの乗算器でk-bitの剰余算ができる(Barrett reductionは(2k, 2k) -> 4k-bitの乗算器を必要とする)
- modが偶数でも動作する(Montgomery乗算はmodが奇数の必要がある)
という二つの性質を持つ。今回は次の6種類のプログラムを実装し、$(8 \times 10^{7})! \bmod 998244353$を計算して速度を比較した。
- 32-bit Montgomery
- 32-bit div2by1
- mod < $2^{30}$という制約を仮定し32-bit div2by1をちょっと改造したもの
- これら三種のアルゴリズムをAVX2で高速化したもの
ベンチマークのコード、およびdiv2by1改造版のpythonコードについては、div2by1 改造版 · GitHub
手元(Ryzen 5 5600X)、およびAtCoderのコードテストでの実行時間は次の通り。
| Local | AtCoder | |
|---|---|---|
| naive | 285ms | 482ms |
| naive(const mod) | 227ms | 298ms |
| montgomery | 202ms | 257ms |
| montgomery AVX | 28ms | 51ms |
| div2by1 | 316ms | 458ms |
| div2by1 AVX | 42ms | 107ms |
| my div2by1 | 292ms | 383ms |
| my div2by1 AVX | 52ms | 83ms |
次のようなことがわかる
- そもそもnaive(=除算命令を使っているはず)が直観よりかなり速い。おそらくIce Lakeから除算が速くなったというやつ(Intel Ice Lakeのプロセッサは整数除算命令がアツい - chroot("/home/hibari"))。でもIce lakeの製造開始は2019年らしい、老人さん？笑
- montgomery + AVX2が速すぎる。4倍速以上になるとは思ってなかった。
- div2by1 あんまり速くない。montgomery速すぎ + AVXで高速化とかするならもうmod 998244353決め打ちで問題なさそう であることを考えると正直出番がなさそうな
