High Performance Correctly Rounded Math Libraries for 32-bit Floating Point Representations PLDI2021読み会@Idein 2021/09/24 藤井裕也

概要 ● 正しく丸められた32bitの数学ライブラリを実装する ○ 正しく丸められた、とは ■ 正しい実数から32bitに丸めた値のこと ○ よく使われているライブラリには誤った出力がある ○ 正確なライブラリもあるが性能が低い ● 著者らの16bit版のライブラリが存在する ○ 32bitにそのまま適用するにはいくつか課題がある ● この論文では16bit版の手法を改良して32bit版を実装する ○ floatとposit32の両方 ■ posit32: 0付近の精度が高いような浮動小数点の規格

先行研究の問題点 ● mini-max approach ○ すべての入力に対する誤差のうち、最大の誤差を最小化する ○ 誤差発生ポイントがいくつかある ■ range reduction ■ output compensation ■ 多項式の評価 ● RLIBM ○ 16bit版の著者らの実装 ○ 上記の誤差発生ポイントを抑えて、結果的に誤差がないようにしている ○ 線形計画問題にしてLP solverで解いている ■ 32bitだと制約の数が多すぎる

論文のキモ ● 16bit版の手法を改良して、LP solverで解けるようにする ● いろんな手段でLPの制約の数を減らしていく ○ counterexample guided polynomial generation ○ piecewise polynominal ○ range reductionの改良

改良点1: counterexample guided polynomial generation ● 入力を全部試さず、サンプリングする ● サンプリングした入力でLPを解いて多項式を生成して、全入 力をテスト ● テストに落ちる入力があればサンプルに加えてLPを解く...を 繰り返し

改良点2: piecewise polynomials ● 入力をsub-domainに分割して、各sub-domainに対して個別 の多項式を生成 ● 入力の範囲とそれに対応する多項式を表で持っておく

改良点3: range reductionの改良 ● 実装する数学関数の特性を使う ○ そのため、関数ごとに別の方法を考えなければならない ● 論文中ではsinの例で説明している ○ 他の例は別論文になっているらしい

sin(𝛑x)の例 ● 32bit = 40億通り以上の入力 ● 32bit floatではsin(𝛑x)は以下の3つの範囲に分けられる ● 一番上のパターンだけで8億通り近くあるので、range reductionを使う

sin(𝛑x)の例 range reduction (1/2) ● sinの三角関数や周期関数の特性を使って考慮すべき入力を 減らす ● ○ ○ J = K + L (Kは整数部分、Lは小数部分) ● [0.5, 1)は[0, 0.5)の鏡像なので、入力は[0,0.5)だけ考えれば 良い ● これでもまだ1億8400万通りくらい入力がある

sin(𝛑x)の例 range reduction (2/2) ● ○ L’ = N/512 + R ● sinpi(N/512)とcospi(N/512)は事前に計算しておく ● ● 入力は0から1/512の範囲まで減らせたが、cosも必要になっ た

sin(𝛑x)の例 多項式の生成 (1/5) 1. 正しいライブラリを使って、各入力に対する正解値を求める 2. range reductionの結果狭くなった入力Rに対する正解値を求 める 3. 入力をsub-domainに分割する 4. 各sub-domainに対する多項式を生成する

sin(𝛑x)の例 多項式の生成 (2/5) ● 正しいライブラリを使って、各入力に対する正解値の範囲を 求める ○ rounding intervalと呼んでいる ● まず32bitでの正解値を求める ● doubleの精度で、丸めると正解値になる上限と下限を求める

sin(𝛑x)の例 多項式の生成 (3/5) ● range reductionの結果の入力Rに対する正解値を求める ● 範囲は狭まったが、sinpi(R)とcospi(R) (R ∈[0,1/512)) が必 要になった ● 32bitで正しい結果になるsinpi(R)とcospi(R)を求める ○ sinpi(R)とcospi(R)の結果が32bitで正しいことは不十分 ○ sinpi(R)とcospi(R)を使った結果の式が32bitで正しい必要がある

sin(𝛑x)の例 多項式の生成 (4/5) ● 入力をsub-domainに分割する ○ 入力を減らしたとはいえ、まだ1億1千万通りもある ○ ビットパターンで入力を分割する ● sub-domainのインデックスのbitはMSBから始めない ○ 事前に範囲を絞ったので、R=0以外の入力は左から6bit分共通している

sin(𝛑x)の例 多項式の生成 (5/5) ● 各sub-domainに対する多項式を生成する 1. 入力をいくつかサンプルしてLP solverで解く 2. 得られた多項式でsub-domainの入力全部をテストする 3. 正しくない答えになった入力があれば、サンプルに加えて1に 戻る

精度の評価 (32bit float)

精度の評価 (posit32) ● posit32向けのライブラリは無いが、精度が良い0付近でも doubleと同じ ● doubleのライブラリと比較

多項式生成にかかった時間 (float) 2.10GHz Intel Xeon Gold 6230Rで実行

多項式生成にかかった時間 (posit32)

実行時間の評価 ● 平均してglibcの1.1倍、intelの1.5倍速い ○ ベクトル化するとintelのライブラリのほうが速くなるが誤差がすごい増えるらし い

まとめ ● 32bitの数学ライブラリを作成した ○ 正確に丸められている ○ 既存のライブラリより高速に動作する ○ posit32向けのライブラリは初