CPU/GPU高速化セミナー浮動小数点から文字列への高速変換の論文を読んでみた / cpugpu acceleration seminar 20230201

by 株式会社フィックスターズ

Slide 1

Slide 1 text

Slide 2

Slide 2 text

Slide 3

Slide 3 text

Slide 4

Slide 4 text

Copyright© Fixstars Group 本講演の位置づけ 4 ● 弊社がサービス展開している、ソフトウェアによる処理速度向上に関連した様々な技術情報を発信しています ● 性能モデルの理論と実践 (発表資料) ○ ソフトウェア性能の理論的根拠となる性能モデルの使用方法を実践例 ● 今回の内容 ○ Ryū のアルゴリズムについて、理論と実装を紹介 ○ GPU移植の説明と、CPU との処理速度を比較 ● こんな方に向いています ○ アルゴリズム・高速計算の実装に興味のある方 ○ システム性能向上に向けた、技術アプローチに興味のある方

Slide 5

Slide 5 text

Copyright© Fixstars Group 発表者紹介冨田明彦ソリューションカンパニー営業企画 2008年に入社。金融、医療業界において、ソフトウェア高速化業務に携わる。その後、新規事業企画、半導体業界の事業を担当し、現職。 5 大野真暉ソリューション第一事業部シニアエンジニア 2019年に新卒で入社。製造業におけるスケジューリング問題のソルバー開発に従事。数理最適化、量子コンピュータ関連のプロジェクトを手掛ける。

Slide 6

Slide 6 text

Slide 7

Slide 7 text

Copyright© Fixstars Group フィックスターズの強みコンピュータの性能を最大限に引き出す、ソフトウェア高速化のエキスパート集団ハードウェアの知見アルゴリズム実装力各産業・研究分野の知見 7 目的の製品に最適なハードウェアを見抜き、その性能をフル活用するソフトウェアを開発します。ハードウェアの特徴と製品要求仕様に合わせて、アルゴリズムを改良して高速化を実現します。開発したい製品に使える技術を見抜き、実際に動作する実装までトータルにサポートします。

Slide 8

Slide 8 text

Copyright© Fixstars Group サービス提供分野 8 半導体自動車産業機器生命科学金融 ● NAND型フラッシュメモリ向けファームウェア開発 ● 次世代AIチップの開発環境基盤 ● 自動運転の高性能化、実用化 ● 次世代パーソナルモビリティの研究開発 ● Smart Factory実現への支援 ● マシンビジョンシステムの高速化 ● ゲノム解析の高速化 ● 医用画像処理の高速化 ● AI画像診断システムの研究開発 ● デリバティブシステムの高速化 ● HFT(アルゴリズムトレード)の高速化

Slide 9

Slide 9 text

Slide 10

Slide 10 text

Copyright© Fixstars Group コンサルティングセットアップ支援、処理の分割や変換等のコンサルクラウド実行環境のご提供クラウド経由での量子コンピュータ利用サービスを提供ソフトウェア高速化・開発支援サービス量子コンピュータを組み合わせてシステムの高速化を実現セミナー・トレーニング量子コンピュータの研究動向や活用事例、実際の利用方法等量子コンピューティング活用支援多様なハードウェアでのソフトウェア高速化サービスに加え、量子コンピュータ活用支援とシステム開発を提供しています。お客様の課題ご支援内容量子コンピューティングが課題の解決に役立つか、確信が持てない量子コンピューティングの検討をどう進めていったら良いかわからない作りたいアプリケーションがあるが、開発が難しい 10

Slide 11

Slide 11 text

Copyright© Fixstars Group 物流倉庫における人員配置最適化の社会実装開発内閣府主催SIP公開シンポジウム2021 11 1 通販向け物流倉庫で誰がどの仕事を担当するか割り振る量子コンピューティングシステムを開発 2 30分以上かかっていた割り当て作業が1秒未満の計算で実現 3 業務負荷の偏りが解消し、全体の作業効率も30％改善する割り当てパターンが計算できた作業者の負荷のばらつきを抑えつつ、全体の作業効率は30%向上物流分野量子コンピューティング活用支援サービス領域ベルメゾンロジスコ様のご協力を得て実施

Slide 12

Slide 12 text

Slide 13

Slide 13 text

Copyright© Fixstars Group Ryū 2018年のPLDIで発表された論文 “Ryū: fast float-to-string conversion” を紹介します。浮動小数点を高速に文字列に変換するアルゴリズムを提案・実装した論文です。固定幅整数の演算のみで文字列に変換し、正確かつ高速なことが特徴です。 Adams, Ulf. "Ryū: fast float-to-string conversion." Proceedings of the 39th ACM SIGPLAN Conference on Programming Language Design and Implementation. 2018. 13

Slide 14

Slide 14 text

Slide 15

Slide 15 text

Slide 16

Slide 16 text

Copyright© Fixstars Group 浮動小数点についての復習：IEEE754（１） IEEE754は浮動小数点数のビット表現（交換形式）を定義している。ビット表現は三つの要素から構成される。浮動小数点の形式ごとに、それぞれの要素のビット幅が決まっている。 1. 符号 sign 2. 指数 exponent 3. 仮数 fraction（mantissaとも） 16 Float64のビット表現 [5]より引用符号のビット幅=1、指数のビット幅=11、仮数のビット幅=52

Slide 17

Slide 17 text

Copyright© Fixstars Group 浮動小数点についての復習：IEEE754（２）浮動小数点のビット表現と数値の対応関係は下式で定義される。 • 𝑠：符号部の値 • 𝑠 = 0なら正の値、𝑠 = 1なら負の値 • 𝑒：指数部の値 • biasは浮動小数点の形式ごとに決まっている値 • Float32なら126で、Float64なら1022 • 𝑚：仮数部のビット列 • 1.mと0.mは2進数で記載した小数 17 浮動小数点の数値 [1]より引用

Slide 18

Slide 18 text

Copyright© Fixstars Group 文字列化のルール浮動小数点をどのような文字列に変換するか？本発表では、次のルール1, 2に従う、正確かつ短い文字列への変換を考える。 a. 情報の保存：パーサ3は文字列化された値から元々の浮動小数点を復元できる b. 文字列は最短：文字列は可能な限り最短にする c. 正確な丸め：aとbのルールを満たす範囲で文字列は浮動小数点の数値に最も近くする 1：How to print floating-point numbers accurately Guy L. Steele, Jr. and Jon L. White. 1990. 2：SteeleとWhiteは他のルールも定めているが、この三つのルールのみ用いられているケースが多い。 3：パーサは文字列化の逆向きの変換つまり、文字列から浮動小数点に変換するアルゴリズムのこと。 18

Slide 19

Slide 19 text

Copyright© Fixstars Group アルゴリズム概要 ~1.2を具体例に~（１） 1.2を例に文字列化アルゴリズムの概要を紹介する。 1.2のFloat64のビット表現は次のようになる。 0011111111110011001100110011001100110011001100110011001100110011 このビット表現を定義通り数値に変換すると次の値になる。 𝑓 =1.1999999999999999555910790149937383830547332763671875 この文字列は正確だがあまりに長い。先ほどのルールの「2. 文字列は最短」に違反する。 19 符号指数仮数浮動小数点の数値 [1]より引用

Slide 20

Slide 20 text

Copyright© Fixstars Group アルゴリズム概要 ~1.2を具体例に~（２）数値𝒇の「近く」で、桁数の少ない数値を探す。「近く」を正確に表現するために、下の数直線のように隣の浮動小数点数との境界を考える。 𝑓− ： 𝑓より一つ小さなFloat64 𝑓+ ： 𝑓より一つ大きなFloat64 𝑓の区間の中で、十進法の小数で表現したときに最も桁数の少ない実数が1.2なので、”1.2”と文字列化する。この文字列化アルゴリズムは先ほど定めた文字列化のルールに従う。 20 𝑓−の区間 𝑓− 𝑓+の区間 𝑓+ 𝑓の区間 𝑓 1.2 𝑓−と𝑓の中点 1.19999999999999984456877655247808434069156646728515625 𝑓+と𝑓の中点 1.20000000000000006661338147750939242541790008544921875

Slide 21

Slide 21 text

Slide 22

Slide 22 text

Copyright© Fixstars Group おことわり正の浮動小数点の文字列化のみを考える。つまり、以下の入力は除外する。 • NaN • Infinity • 負の浮動小数点 • 負の値を文字列化するときは、その正の値の文字列の先頭に「-」を付けると良い • +0.0 • アルゴリズムの都合上除外するまた、指数部の値が変わる特殊なケースや丸めなどの話はスキップする。（大事ではあるがアルゴリズムの概要を把握する分には不要） 22

Slide 23

Slide 23 text

Copyright© Fixstars Group 文字列化アルゴリズム：理論面 1. 浮動小数点を数値にデコードする 2. 隣の浮動小数点数との中点を計算する ○ 小さな浮動小数点数との中点、元の浮動小数点数、大きな浮動小数点数との中点の三つ組を 𝑢, 𝑣, 𝑤 × 2𝑒2とする（𝑢, 𝑣, 𝑤は正の整数、 𝑒2 は整数） 3. 2の冪を10の冪となるよう計算する ○ 三つ組を 𝑎, 𝑏, 𝑐 × 10𝑒10の形に変形する（𝑎, 𝑏, 𝑐は正の整数、 𝑒10 は整数） 4. 主要部の桁数が最も小さく正しく丸められた整数を計算する ○ 主要部とは、整数をm × 10𝑒と書いた時のmのこととする ○ 𝑎より大きく𝑐より小さい整数の中で主要部の桁数が小さい整数𝑑𝑜 , 𝑒𝑜 を計算する 5. 整数を文字列としてプリントする ○ 𝑑𝑜 × 10𝑒𝑜+𝑒10を文字列にする 23

Slide 24

Slide 24 text

Slide 25

Slide 25 text

Slide 26

Slide 26 text

Copyright© Fixstars Group STEP2：隣の浮動小数点数との中点を計算する • 概要 • 一つ小さな浮動小数点数𝑓−との中点を計算する • 一つ大きな浮動小数点数𝑓+との中点を計算する • 隣の浮動小数点数はIEEE754の式とにらめっこすれば計算できる • 入力 • 𝑒, 𝑚 • 𝑒𝑓 , 𝑚𝑓 • 出力 • 𝑒2 • 𝑢, 𝑣, 𝑤 26 𝑢 ⋅ 2𝑒2 𝑓+ 𝑓 𝑤 ⋅ 2𝑒2 𝑣 ⋅ 2𝑒2 𝑓− 三つ組の計算式 𝑒2 = 𝑒𝑓 − 2 𝑢 = 4𝑚𝑓 − ቊ 1 if 𝑚 = 0 … 0 2 and 𝑒 > 0 … 01 2 2 otherwise 𝑣 = 4𝑚𝑓 𝑤 = 4𝑚𝑓 + 2 𝑚𝑓 × 2𝑒𝑓 = 𝑣 × 2𝑒2 が成立する。

Slide 27

Slide 27 text

Copyright© Fixstars Group STEP2：具体例 1.2 27 𝑢 ⋅ 2𝑒2 𝑓+ 𝑓 = 𝑚𝑓 × 2𝑒𝑓 𝑤 ⋅ 2𝑒2 𝑣 ⋅ 2𝑒2 𝑓− 𝑒2 = 𝑒𝑓 − 2 = −54 𝑣 = 4 × 𝑚𝑓 = 21617278211378380 𝑢 = 𝑣 − 2 = 21617278211378378 𝑤 = 𝑣 + 2 = 21617278211378382 2. 隣の浮動小数点との中点を計算する

Slide 28

Slide 28 text

Copyright© Fixstars Group STEP3：2の冪を10の冪となるよう計算する • 概要 • 10進法の数値として文字列化する際に2𝑛は邪魔なので展開する • 入力 • 𝑒2 • 𝑢, 𝑣, 𝑤 • 出力 • 𝑒10 • 𝑎, 𝑏, 𝑐 28 𝑢 × 2𝑒2 𝑓 = 𝑣 × 2𝑒2 𝑤 × 2𝑒2 𝑎 × 10𝑒10 𝑏 × 10𝑒10 𝑐 × 10𝑒10 三つ組（10進数）の計算式 𝑒2 ≥ 0 ⇒ 𝑒10 = 0, 𝑎, 𝑏, 𝑐 = 𝑢, 𝑣, 𝑤 ⋅ 2𝑒2 𝑒2 < 0 ⇒ 𝑒10 = 𝑒2 , 𝑎, 𝑏, 𝑐 = 𝑢, 𝑣, 𝑤 ⋅ 5−𝑒2 𝑎, 𝑏, 𝑐 × 10𝑒10 = 𝑢, 𝑣, 𝑤 × 2𝑒2 が成立する。

Slide 29

Slide 29 text

Copyright© Fixstars Group STEP3：具体例 1.2 29 𝑢 × 2𝑒2 𝑣 × 2𝑒2 𝑤 × 2𝑒2 𝑎 × 10𝑒10 𝑏 × 10𝑒10 𝑐 × 10𝑒10 𝑒2 ≥ 0 ⇒ 𝑒10 = 0, 𝑎, 𝑏, 𝑐 = 𝑢, 𝑣, 𝑤 ⋅ 2𝑒2 𝑒2 < 0 ⇒ 𝑒10 = 𝑒2 , 𝑎, 𝑏, 𝑐 = 𝑢, 𝑣, 𝑤 ⋅ 5−𝑒2 𝑒10 = 𝑒2 = −54 𝑎 = 𝑢 × 554 = 1199999999999999844568776552478084340691566467285156250 𝑏 = 𝑣 × 554 = 1199999999999999955591079014993738383054733276367187500 𝑐 = 𝑤 × 554 = 1200000000000000066613381477509392425417900085449218750 3. 2の冪を10の冪となるよう計算する

Slide 30

Slide 30 text

Copyright© Fixstars Group STEP4：主要部の桁数が最も小さく正しく丸められた整数を計算する • 概要 • compute_shortestで、𝑎より大きく𝑐より小さい整数の中で、主要部の桁数が最も小さい整数を求める • 数学的に定式化すると、次の不等式を満たす整数𝑒𝑜 , 𝑑𝑜 を求めることに対応する • 𝒂 < 𝒅𝒐 × 𝟏𝟎𝒆𝒐 < 𝒄 • ただし𝑑𝑜 の桁数は最小とする • つまり𝑒𝑜 の値は最大 • 入力 • 𝑎, 𝑏, 𝑐 • 丸めのオプション • 出力 • 𝑒𝑜 , 𝑑𝑜 30 𝑎と𝑐が等しくなるまで 10で除算する compute_shortest [1]より引用

Slide 31

Slide 31 text

Copyright© Fixstars Group STEP4：（資料）正しいcompute_shortest 右の関数が真のcompute_shortestである。前ページのcompute_shortestには一部嘘があった。嘘とは、文字列化のルールの一つに違反していること。 (3) 正確な丸め：1と2のルールを満たす文字列で、文字列は浮動小数点の数値に最も近くする隣の浮動小数点数との中点である𝑎と𝑐のみならず、𝑏も一緒に計算することでこのルールも満たしたアルゴリズムになる。 31 def compute_shortest(a, b, c, accept_smaller, accept_larger, break_tie_down): i=0 a[0]=a, b[0]=b, c[0]=accept_larger ? c : c-1 all_a_zero[0]=True, all_b_zero[0]=True digit[0]=0 while floor(a[i]/10) < floor(c[i]/10): all_a_zero[i+1]=all_a_zero[i] and a[i]%10==0 all_b_zero[i+1]=all_b_zero[i] and digit[i]==0 a[i+1]=floor(a[i]/10), b[i+1]=floor(b[i]/10), c[i+1]=floor(c[i]/10) digit[i+1]=b[i]%10 i=i+1 if accept_smaller and all_a_zero[i]: while a[i]%10==0: all_b_zero[i+1]=all_b_zero[i] and digit[i]==0 a[i+1]=a[i]/10, b[i+1]=floor(b[i]/10), c[i+1]=floor(c[i]/10) digit[i+1]=b[i]%10 i=i+1 is_tie=(digit[i]==5) and all_b_zero[i] wand_round_down=(digit[i]<5) or (is_tie and break_tie_down) round_down=(want_round_down and (a[i]!=b[i] or all_a_zero[i])) or (b[i]+1>c[i]) return (d_o, e_o) = round_down ? (b_i, i) : (b_i+1, i)

Slide 32

Slide 32 text

Slide 33

Slide 33 text

Slide 34

Slide 34 text

Slide 35

Slide 35 text

Copyright© Fixstars Group 文字列化アルゴリズム：理論面（再掲） 1. 浮動小数点を数値にデコードする 2. 隣の浮動小数点数との中点を計算する ○ 小さな浮動小数点数との中点、元の浮動小数点数、大きな浮動小数点数との中点の三つ組を 𝑢, 𝑣, 𝑤 × 2𝑒2とする（𝑢, 𝑣, 𝑤は正の整数、 𝑒2 は整数） 3. 2の冪を10の冪となるよう計算する ○ 三つ組を 𝑎, 𝑏, 𝑐 × 10𝑒10の形に変形する（𝑎, 𝑏, 𝑐は正の整数、 𝑒10 は整数） 4. 主要部の桁数が最も小さく正しく丸められた整数を計算する ○ 主要部とは、整数をm × 10𝑒と書いた時のmのこととする ○ 𝑎より大きく𝑐より小さい整数の中で主要部の桁数が小さい整数𝑑𝑜 , 𝑒𝑜 を計算する 5. 整数を文字列としてプリントする ○ 𝑑𝑜 × 10𝑒𝑜+𝑒10を文字列にする 35

Slide 36

Slide 36 text

Copyright© Fixstars Group 文字列化アルゴリズム：具体例（まとめ） 36 S T E P 1 0 s 1023 e 900719925474099 m 𝑓 = −1 0 × 252 + 900719925474099 × 21023 −1023 −52 𝑚𝑓 = 5404319552844595 𝑒𝑓 = −52 𝑢 ⋅ 2𝑒2 𝑓+ 𝑓 𝑤 ⋅ 2𝑒2 𝑣 ⋅ 2𝑒2 𝑓− 𝑒2 = 𝑒𝑓 − 2 = −54 𝑣 = 4 × 𝑚𝑓 = 21617278211378380 𝑢 = 𝑣 − 2 = 21617278211378378 𝑤 = 𝑣 + 2 = 21617278211378382 S T E P 2

Slide 37

Slide 37 text

Copyright© Fixstars Group 文字列化アルゴリズム：具体例（まとめ） 37 𝑢 × 2𝑒2 𝑣 × 2𝑒2 𝑤 × 2𝑒2 𝑎 × 10𝑒10 𝑏 × 10𝑒10 𝑐 × 10𝑒10 𝑒2 ≥ 0 ⇒ 𝑒10 = 0, 𝑎, 𝑏, 𝑐 = 𝑢, 𝑣, 𝑤 ⋅ 2𝑒2 𝑒2 < 0 ⇒ 𝑒10 = 𝑒2 , 𝑎, 𝑏, 𝑐 = 𝑢, 𝑣, 𝑤 ⋅ 5−𝑒2 𝑒10 = 𝑒2 = −54 𝑎 = 𝑢 × 554 = 1199999999999999844568776552478084340691566467285156250 𝑏 = 𝑣 × 554 = 1199999999999999955591079014993738383054733276367187500 𝑐 = 𝑤 × 554 = 1200000000000000066613381477509392425417900085449218750 S T E P 3 S T E P 4 S T E P 5 ⋮ 𝑎 / 1052 = 119 < 𝑐 / 1052 = 120 𝑎 / 1053 = 11 < 𝑐 / 1053 = 12 𝑎 / 1054 = 1 == 𝑐 / 1054 = 1 ⇒ 𝑑𝑜 = 12, 𝑒𝑜 = 53 𝑑𝑜 × 10𝑒𝑜+𝑒10を文字列にする 𝑑𝑜 = 12 𝑒𝑜 + 𝑒10 = −1 ⇒”1.2”

Slide 38

Slide 38 text

Copyright© Fixstars Group 文字列化アルゴリズム：実装面先ほどのアルゴリズムのボトルネックは任意精度の演算を行う部分である。任意精度の演算を行うのはSTEP3とSTEP4である。 • STEP3：2の冪を処理して浮動小数点を「整数×10の冪」と表す • この整数が大きな整数となる • STEP4：二つの大きな整数が等しくなるまで10で除算する 38 1. 浮動小数点を数値にデコードする 2. 隣の浮動小数点数との中点を計算する ○ 小さな浮動小数点数との中点、元の浮動小数点数、大きな浮動小数点数との中点の三つ組を 𝑢, 𝑣, 𝑤 × 2𝑒2とする（𝑢, 𝑣, 𝑤は正の整数、 𝑒2 は整数） 3. 2の冪を10の冪となるよう計算する ○ 三つ組を 𝑎, 𝑏, 𝑐 × 10𝑒10の形に変形する（𝑎, 𝑏, 𝑐は正の整数、 𝑒10 は整数） 4. 主要部の桁数が最も小さく正しく丸められた整数を計算する ○ 主要部とは、整数をm × 10𝑒と書いた時のmのこととする ○ 𝑎より大きく𝑐より小さい整数の中で主要部の桁数が小さい整数𝑑𝑜 , 𝑒𝑜 を計算する 5. 整数を文字列としてプリントする ○ 𝑑𝑜 × 10𝑒𝑜+𝑒10を文字列にする 𝑎, 𝑏, 𝑐 がとても大きな整数

Slide 39

Slide 39 text

Copyright© Fixstars Group 具体例：Float32の𝒂, 𝒃, 𝒄の値 Float32 符号ビット：1、指数ビット：8、仮数ビット：23 bias：127 STEP1 𝑒𝑓 : 0 ~ 224 𝑚𝑓 : −149~104 STEP2 𝑢, 𝑣, 𝑤: 0~ 226 𝑒2 : −151~102 STEP3 𝑎, 𝑏, 𝑐: max 226 × 2102, 226 × 5151 ≒ 10131 39 STEP1 STEP2 STEP3 非常に大きい！

Slide 40

Slide 40 text

Copyright© Fixstars Group 具体例：Float64の𝒂, 𝒃, 𝒄の値 Float32 符号ビット：1、指数ビット：11、仮数ビット：52 bias：1023 STEP1 𝑒𝑓 : 0 ~ 253 𝑚𝑓 : −1074~971 STEP2 𝑢, 𝑣, 𝑤: 0~ 255 𝑒2 : −1076~969 STEP3 𝑎, 𝑏, 𝑐: max 255 × 2969, 255 × 51076 ≒ 10768 40 STEP1 STEP2 STEP3 非常に大きい！

Slide 41

Slide 41 text

Copyright© Fixstars Group 大きな整数の回避（１） STEP3とSTEP4を工夫できないか？ • STEP3：2の冪を処理して浮動小数点数を「整数×10の冪」とあらわす • この整数が大きな整数となる • STEP4：二つの大きな整数が等しくなるまで10で除算する STEP3で「大きな整数」があらわれ、STEP4で「大きな整数を10で除算した小さな整数」にする。 STEP3とSTEP4を繋げることで、「大きな整数」を計算することなく、直接「大きな整数を10で除算した小さな整数」を計算することを考える。 41

Slide 42

Slide 42 text

Copyright© Fixstars Group 大きな整数の回避（２）今後、𝒆𝟐 ≥ 𝟎の場合に限定する。 STEP3とSTEP4で行う計算の概略を右に示す。数学的な考察から、STEP4のwhileループは少なくとも次の𝑖 = 𝑞まではまわることが分かる。 𝑞 = max(0, 𝑒2 × log10 2 − 1) 𝑎を10で𝑞回除算した値𝑎𝑞 は𝑢の100倍未満のため、小さな整数である： 𝑎𝑞 = 𝑎 / 10𝑞 ≒ 𝑢 × 2𝑒2 / 10𝑞 < 100𝑢 よって、whileループが𝑞回まわった時点の変数の値は小さい。 ⇒STEP3の段階で、STEP4のwhileループが𝒊 = 𝒒までまわった時の変数の値を計算できれば、大きな整数の出現を回避できる！以降、 𝑎を計算せずに、直接𝑎𝑞 を計算する方法について考察する。 42 STEP3 𝑎 = 𝑢 × 2𝑒2, 𝑐 = 𝑤 × 2𝑒2 STEP4 while 𝑎 / 10𝑖 < 𝑐 − 1 / 10𝑖: 𝑖 = 𝑖 + 1 2𝑒2 ≒ 10𝑞 となる𝑞を選んでいる

Slide 43

Slide 43 text

Copyright© Fixstars Group 𝒂を10で𝒒回除算した値の計算 compute_shortestに登場する変数は色々あるが、 𝑎を10で𝑞回除算した値𝑎𝑞 の計算が最も困難なので、 𝑎𝑞 の計算方法のみ解説する。 𝑎𝑞 を計算しやすい形に式変形する。大きな値に赤枠を付ける。非整数になりうる値に青枠を付ける。 𝑎𝑞 = 𝑢 × 2𝑒2 / 10𝑞 = 𝑢 × 2𝑒2−𝑞 / 5𝑞 （アイディア：大きな整数2𝑘で5𝑞の除算をなんとかする） = 𝑢 × 2𝑒2−𝑞−𝑘 × 2𝑘 / 5𝑞 （アイディア： 𝑘が十分大きいと次の式が成立する） = 𝑢 × 2𝑒2−𝑞−𝑘 × 2𝑘 / 5𝑞 + 1 43 式変形☆ 次スライド説明する

Slide 44

Slide 44 text

Slide 45

Slide 45 text

Copyright© Fixstars Group 変形した数式について ☆で式変形した下式は以下の点で優れている。実行時変数は𝑢, 𝑒2 , 𝑞, 𝑘である。 • 𝑎𝑞 は三つの項の積である • 赤枠（ 2𝑒2−𝑞−𝑘 ）部分は𝑒2 − 𝑞 − 𝑘 < 0のため除算になるが、シフト演算で効率的に実装できる • 困難な5𝑞による除算が 2𝑘 / 5𝑞 に吸収されている • Float32やFloat64の場合、 𝑞の値はそれほど大きくない • 𝑞ごとに 2𝑘 / 5𝑞 の値をあらかじめ計算してテーブル化すれば、実行時には5𝑞による除算がなくなる赤枠と緑枠の項が大きな値とならないような𝑘が存在すれば、効率的に𝑎𝑞 を計算できる。 ⇒2𝑘 ≒ 5𝑞となる𝑘 = log2 5𝑞付近で、下界より大きな𝑘 は存在するか？ 45 𝑎𝑞 = 𝑢 × 2𝑒2−𝑞−𝑘 × 2𝑘 / 5𝑞 + 1

Slide 46

Slide 46 text

Copyright© Fixstars Group 丁度良い𝒌の探索 • Float64で𝑘の下界を計算した（右図） • 横軸は𝑒2 • 𝑘の下界を青色でプロット • 2𝑘 / 5𝑞 の値をテーブルとして保持するのに必要な要素数の最小値を橙色でプロット • おおよそ 𝑘 − log2 5𝑞 に対応する • 橙色はプラトーになっている！ • 最大値を𝐵0 とする • 𝑘 = 𝐵0 + log2 5𝑞 とすれば、𝑘は下界よりも大きな値となる • さらに、実際に計算すると分かるが、この𝑘は丁度良い𝑘である！（次の二ページで計算する） • 2𝑘と5𝑞は大きな整数だが、 𝑞ごとに緑枠の値をあらかじめ計算しておけば、実行時はテーブルを参照するだけで良い 46 𝑘の下界 [1]より引用

Slide 47

Slide 47 text

Copyright© Fixstars Group Float32の場合 • 変数の整理 • 𝑒2 : 0 ∼ 102, 𝑢: 0 ∼ 226 • 𝑞 = max(0, 𝑒2 × log10 2 − 1) • 𝑘 = 𝐵0 + log2 5𝑞 • 𝐵0 = 60 • 𝑞 + 𝑘 − 𝑒2 ≒ 𝐵0 + 𝑞 × 1 + log2 5 − 𝑒2 ≒ 𝐵0 + 𝑒2 × log10 2 − 1 × log2 10 − 𝑒2 = 𝐵0 + 𝑒2 − log2 10 − 𝑒2 ≒ 𝐵0 = 60 • 2𝑘 / 5𝑞 ≒ 2𝐵0+log2 5𝑞 / 5𝑞 ≒ 2𝐵0 = 260 • 赤枠：60ビット右シフト、緑枠：uint64_tの範囲内 • u（uint32_t）×緑枠（uint64_t）の結果を60ビット右シフトすると𝒂𝒒 （uint32_t）が得られる • uint128_tの範囲で計算できる 47

Slide 48

Slide 48 text

Copyright© Fixstars Group Float64の場合 • 変数の整理 • 𝑒2 : 0 ∼ 969, 𝑢: 0 ∼ 255 • 𝑞 = max(0, 𝑒2 × log10 2 − 1) • 𝑘 = 𝐵0 + log2 5𝑞 • 𝐵0 = 124 • 𝑞 + 𝑘 − 𝑒2 ≒ 𝐵0 + 𝑞 × 1 + log2 5 − 𝑒2 ≒ 𝐵0 + 𝑒2 × log10 2 − 1 × log2 10 − 𝑒2 = 𝐵0 + 𝑒2 − log2 10 − 𝑒2 ≒ 𝐵0 = 124 • 2𝑘 / 5𝑞 ≒ 2𝐵0+log2 5𝑞 / 5𝑞 ≒ 2𝐵0 = 2124 • 赤枠：124ビット右シフト、緑枠：uint128_tの範囲内 • u（uint64_t）×緑枠（uint128_t）の結果を124ビット右シフトすると𝒂𝒒 （uint64_t）が得られる • uint256_tの範囲で計算できる 48

Slide 49

Slide 49 text

Copyright© Fixstars Group 大きな整数の回避：まとめ • STEP3の計算を工夫することで、STEP4のwhileループで𝑖 = 𝑞の時点までスキップし、大きな整数の出現を回避できる • 工夫した計算はFloat32の場合、uint128_tで行える • Float64の場合、uint256_tで行える • 𝑞ごとに下式の緑枠の値を予め計算しておくことで、文字列化の際はテーブル参照するだけで良い 49

Slide 50

Slide 50 text

Copyright© Fixstars Group 文字列化アルゴリズム：実装面（まとめ） 1. 浮動小数点数を数値にデコードする 2. 隣の浮動小数点数との中点を計算する ○ 小さな浮動小数点数との中点、元の浮動小数点数、大きな浮動小数点数との中点の三つ組を 𝑢, 𝑣, 𝑤 × 2𝑒2とする（𝑢, 𝑣, 𝑤は正の整数、 𝑒2 は整数） 3. 2の冪を10の冪となるよう計算する ○ 三つ組を 𝑎, 𝑏, 𝑐 × 10𝑒10の形に変形する（𝑎, 𝑏, 𝑐は正の整数、 𝑒10 は整数） ○ 三つ組を10でq回除算した値(𝑎𝑞 , 𝑏𝑞 , 𝑐𝑞 )を計算する（𝑒10 も計算する） 4. 主要部の桁数が最も小さく正しく丸められた整数を計算する ○ 主要部とは、整数をm × 10𝑒と書いた時のmのこととする ○ 𝑎𝑞 より大きく𝑐𝑞 より小さい整数の中で主要部の桁数が小さい整数𝑑𝑜 , 𝑒𝑜 を計算する 5. 整数を文字列としてプリントする ○ 𝑑𝑜 × 10𝑒𝑜+𝑞+𝑒10を文字列にする 50

Slide 51

Slide 51 text

Copyright© Fixstars Group 性能評価： Ryū 浮動小数点を乱数で生成し、文字列化に要する時間を比較した。 • 横軸：浮動小数点のビット表現を符号なし整数とみなした時の値 • 縦軸：文字列化に要する時間[ns] 横軸の値がFloat32の場合は231の時、Float64の場合は 263 の時、浮動小数点の符号が入れ替わる。 c, dのグラフが線対称なのはこのためである。いずれの場合もRyūが最速だった。 51 浮動小数点を文字列化するのに要する時間の比較 [1]より引用上段：C、下段：Java 左側：Float32、右側：Float64 実験環境 CPU：Intel(R) Core(TM) i7-4770K at 3.50GHz OS：Ubuntu Linux 17.10. C：clang 3.9 with the -O2 flag Java：OpenJDK 1.8.0_131.

Slide 52

Slide 52 text

Copyright© Fixstars Group GPU移植 RyūをGPUに移植し、大量の浮動小数点を文字列化するタスクの高速化を行った。浮動小数点の配列と区切り文字列を受け取り、浮動小数点の文字列を区切り文字でつなげた文字列を出力する関数を実装した。 void toString(const std::vector& v, const std::string& sep, std::string* s); T = float or double 例えばv={1.2, 3.14, -0.7}, sep=”, ”のとき、”1.2, 3.14, -0.7”と文字列化する。 52

Slide 53

Slide 53 text

Copyright© Fixstars Group GPU移植：実装（１） GPUでのtoStringを３つの関数に分けて実装した。 1. f2dKernel, d2dKernel ○ 浮動小数点を正しく丸めた10進数表現を計算する ■ 文字列化アルゴリズムのSTEP1～STEP4に対応する ○ 文字列化したときの文字数も計算する 2. thrust::exclusive_scan ○ prefix sumで、文字列化した浮動小数点を書き込むインデックスを計算する 3. d2sKernel ○ 10進数表現を文字列化し、前ステップで計算したインデックスに文字列を書き込む ■ 文字列化アルゴリズムのSTEP5に対応する 53 float to decimal double to decimal という意味

Slide 54

Slide 54 text

Copyright© Fixstars Group GPU移植：実装（２） 54 1.2 3.14 -0.7 0.006 … v 12 314 -7 6 … f2dKernel, d2dKernel GPUの各スレッドが配列の各要素を処理する -1 -2 -1 -3 … d e 3 4 4 5 … size 0 3 7 11 16 … index thrust::exclusive_scan size配列のprefix sumを計算する 1 . 2 3 . 1 4 - 0 . 7 0 . 0 0 6 … s d2sKernel GPUの各スレッドが配列の各要素を文字列化し、 indexで指定されたインデックスに書き込む区切り文字列が空文字列の場合

Slide 55

Slide 55 text

Copyright© Fixstars Group GPU移植：性能評価（１） 3つのアルゴリズムを比較した。 1. Ryū（CPU） 2. Ryū（GPU） 3. yyjson（CPU） ○ Cで実装された、高速なJSONライブラリ ○ 浮動小数点の文字列化はSchubfachベースのアルゴリズムで行う ○ ※浮動小数点配列を文字列化した結果はRyūと異なる ○ ※フェアな比較対象ではないが、参考として比較した Ryū（CPU）とyyjson（CPU）はCPUのシングルスレッドで動作するプログラムである。 55

Slide 56

Slide 56 text

Copyright© Fixstars Group GPU移植：性能評価（２）実験環境 ○ CPU：Intel(R) Xeon(R) CPU E5-2698 v4 @ 2.20GHz ○ GPU：Tesla V100 計測区間 1. Ryū（CPU）：配列を文字列化する時間 2. Ryū（GPU）：GPU上の配列をGPU上で文字列化する時間 3. yyjson（CPU）：文字列のマロック＋JSONオブジェクトを文字列化する時間 ○ yyjson_mut_writeの時間 56

Slide 57

Slide 57 text

Copyright© Fixstars Group GPU移植：性能評価（３） • Ryū（CPU）は浮動小数点を一要素当たり、floatは43[ns]で、doubleは60[ns]で文字列化した • おおよそ論文通りのオーダーだった • Ryū（GPU）はCPU版よりfloatで200倍ほど、doubleで100倍ほど高速となった • 実行時間の8割はd2sKernelで、これを高速化すればもっと速くなる • yyjson（CPU）はRyū（CPU）よりも高速だった • 本実験でも、SchubfachはRyūよりも高速だった 57 配列の要素数 Ryū（CPU） Ryū（GPU） 1,000,000 46.2 0.513 10,000,000 433 2.42 100,000,000 4,335 22.2 配列の要素数 Ryū（CPU） Ryū（GPU） yyjson（CPU） 1,000,000 60.3 0.795 53.7 10,000,000 604 5.60 533 100,000,000 6,037 54.1 5,280 floatの配列を文字列化するのに要する時間 [ms] yyjson（CPU）は浮動小数点を全てdoubleで扱うため除外 doubleの配列を文字列化するのに要する時間 [ms]

Slide 58

Slide 58 text

Slide 59

Slide 59 text

Copyright© Fixstars Group 最新アルゴリズムの紹介：比較 59 初出時期著者アルゴリズム名前の由来 Schubfach OpenJDKのcore- libs-devのパッチ 2018年4月, 11 月 Raffaello Giulietti Ryūのほぼ上位互換 Ryūのcompute_shortest をループ無しで行う鳩ノ巣原理のドイツ語の Schubfachprinzipに由来する Ryū PLDI 2018年6月 Ulf Adams 本発表で解説した Grisuの別名のDragonを日本語に訳した Grisu-Exact GitHub 2020年 jk-jeon Ryū + Grisu Dragonbox GitHub 2020年 jk-jeon Schubfach + Grisu Grisu⇒Dragon Schubfach⇒box (鳩ノ巣原理＝箱入れ原理) ※Grisu-ExactとDragonboxのアルゴリズムの「+ Grisu」について RyūとSchubfachは三つ組で計算する。 Float64の場合、64bit×128bitを三回も計算する必要があり、高コストである。 Grisuのアイディアを利用して、三つ組の計算を省略し、一度の64bit×128bitの計算で済ませるようにした。

Slide 60

Slide 60 text

Copyright© Fixstars Group 最新アルゴリズムの紹介：性能評価 • 文字列化に要する時間を4つのアルゴリズムで比較した ([6]より引用) • 上段：Float32、下段：Float64 • 左側：最短文字列がN桁の浮動小数点を文字列化するのに要する時間（Nが横軸） • 右側：浮動小数点をビットとして一様に生成し、文字列化に要する時間をプロット • 結果 • Dragonboxが最速だった • 最短文字列の桁数が小さいほど、Ryūは他のアルゴリズムよりも遅かった • compute_shortest でのナイーブなループが遅い 60 実験環境 CPU：Intel (R) Core™ i7- 7700HQ CPU @2.80GHz C++：clang-cl compiler shipped with Visual Studio 2022 17.0.4 凡例を拡大

Slide 61

Slide 61 text

Copyright© Fixstars Group 参考資料 • [1] https://dl.acm.org/doi/10.1145/3192366.3192369 • Ryūの論文 • [2] https://github.com/ulfjack/ryu • Ryū著者による実装 • [3] https://github.com/abolz/Drachennest • 浮動小数点をminimal length decimalに変換するアルゴリズムの実行時間を比較する • 文字列としてプリントするのに要する時間を除いて計測し、比較している • [4] https://kubo39.hatenablog.com/entry/2022/06/09/浮動小数点数の10進表記についてまとめてみた • 各種文字列化アルゴリズムがどのライブラリで使用されているか調べている • [5] https://ja.wikipedia.org/wiki/IEEE_754 • Wikipedia IEEE754 • [6] https://github.com/jk-jeon/dragonbox/blob/master/other_files/Dragonbox.pdf • Dragonboxのアルゴリズムを説明したpdf • [7] https://github.com/ibireme/yyjson • yyjson 61