コグノスケ

三角関数のfloat版

目次: C言語とlibc

標準Cライブラリにはdoubleを返す三角関数（sin(), cos(), tan()）とfloatを返す三角関数（sinf(), cosf(), tanf()）が定義されています。

標準Cライブラリの一つの実装であるmuslのコードを見ると、sinf, cosf, tanfの計算にdouble演算を内部で使っています。これは基になったFreeBSDのlibmと同じ実装です。PCはfloatでもdoubleでも関係なく速いんですが、doubleをハードで扱えない貧相なプロセッサには優しくない作りです。

もう一つの実装であるnewlibのコードを見ると、double版のsin, cosこそFreeBSDの実装と同じですが、float版のsinf, cosfはfloatを使ったコードが独自に追加されていて、貧相なプロセッサにも優しい作りになっています。組み込みにやたら使われる実績は伊達じゃないですね。

じゃあmuslにnewlibのfloat版のsinf, cosfを移植すれば、doubleが苦手なプロセッサでも速くなるのでは？と思いました。

テストを先に書こう

コードを触る前に、それぞれの実装の素性を調べておこうと思います。テスト方法は、

• 期待値: glibcのsin, cos, tan（double版）をfloatに変換した結果
• 判定方法: sinf, cosfは1の誤差を許す、tanfは3の誤差を許す
• 比較範囲: floatを32bitとして、全値域（＝約43億パターン）

どうしてtanfだけ判定が甘いかというと、正しい値がわからなかったからです。なぜかglibcの実装も誤差1に収まっていません。どういうことなの……。

テストのコードはGitHubに置き（リンク）ました。特に難しい点はありませんが、muslとnewlibから三角関数を拝借するところは、やや面倒かもしれません。

現在のプロセッサは超速いし、問題の性質上マルチスレッド化も簡単ですから、32並列くらいで頑張れば32bit全域を調査しても3分もかかりません。楽勝ですね〜。

最初にテストして良かった

テスト結果は、当然、全て一致かと思いきや、そんなことはなかった。最初にテストしておいて良かったですね。

良い方から言うとmuslは内部でdoubleで演算しているからか、結果もパーフェクトでした。

一方のnewlibはcosfだけ変な値を返します。32bit全域を試してわずか6パターンです。

cos,cosf_newlib: NG : x:3fc90fe0 f:1.570797 d:1, exp:b52bbbd3 -0.000001, res:b52bbbd0 -0.000001
cos,cosf_newlib: NG : x:3fc90fe1 f:1.570797 d:1, exp:b54bbbd3 -0.000001, res:b54bbbd0 -0.000001
cos,cosf_newlib: NG : x:3fc90fe2 f:1.570797 d:1, exp:b56bbbd3 -0.000001, res:b56bbbd0 -0.000001

cos,cosf_newlib: NG : x:bfc90fe0 f:-1.570797 d:1, exp:b52bbbd3 -0.000001, res:b52bbbd0 -0.000001
cos,cosf_newlib: NG : x:bfc90fe1 f:-1.570797 d:1, exp:b54bbbd3 -0.000001, res:b54bbbd0 -0.000001
cos,cosf_newlib: NG : x:bfc90fe2 f:-1.570797 d:1, exp:b56bbbd3 -0.000001, res:b56bbbd0 -0.000001

正負を考慮（浮動小数点は最上位ビットが符号を示すビット）すると、実質3パターンで変な値が返ることがわかります。

• 正: 0x3fc90fe0, 0x3fc90fe1, 0x3fc90fe2
• 負: 0xbfc90fe0, 0xbfc90fe1, 0xbfc90fe2

誤差は3でした。ほぼ合ってます、おしい。誤差が出ることも不思議ですが、sinfは合っていてcosfだけ値がズレるのも不思議です。