✊🏻 🏂🏼 👨‍⚕️ 正弦関数の例を使用した浮動小数点数の正確で高速な計算。パート2：libm 💪🏼 😥 🔞

フローティングポイントの操作に関する一連の記事を続けます。では最初の記事私は小さな数学的な導入を与え、さまざまな「落とし穴」とプログラムの例の正弦を計算するための最も簡単かつ明白な方法を示しました。今日の記事はスタイルが少し異なります。ここでは練習はありませんが、数学を深く掘り下げて、聖なるもの、つまり標準のライブラリコードについて説明します。また、最初の記事の最後にある質問に答えます。じゃあ行こう。

もう一度いくつかの数学

明らかに、テイラー系列の絶対値が速く減少するほど、必要な精度を達成するために必要な項は少なくなります。そのため、結果はより正確になるようです（これについては以下で詳しく説明します）。比較のために、たとえば、テイラーシリーズの7度の項を取ります（

\frac{x^{7}}{7!}

$\frac{x^7}{7!}$ ）で

x = 1.0

$x=1.0$ そして

x = 0.1

$x=0.1$ ..。式の値は次のようになります

1.198 \cdot 10^{- 4}

$1.198\cdot 10^{-4}$ そして

1.198 \cdot 10^{- 11}

$1.198\cdot 10^{-11}$ それぞれ。大きな違いですね。それでは、正弦関数の計算間隔の上限を減らす方法を見つけてみましょう。

与えられた値の周りのシリーズ展開

この方法を理解するには、研究所の1年目に戻り、テイラーシリーズ（wiki）の定義を覚えておく必要があります。一言で言えば：ある時点で関数とその派生物を知っていると、テイラーシリーズに展開することで、この点の近くにある関数の値を見つけることができます。正弦関数の場合、これは次のことを意味します

\sin (x_{0} + Δ x) \approx \sin (x_{0}) + \sin^{'} (x_{0}) \frac{Δ x}{1!} + \sin^{″} (x_{0}) \frac{Δ x^{2}}{2!} + \sin^{‴} (x_{0}) \frac{Δ x^{3}}{3!} + \dots

$\sin(x_0+\Delta x) \approx \sin(x_0)+\sin'(x_0)\frac{\Delta x}{1!}+\sin''(x_0)\frac{\Delta x^2}{2!}+\sin'''(x_0)\frac{\Delta x^3}{3!}+\cdots$

このアプローチは、実用的な観点から私たちに何をもたらしますか？からの間隔があると想像してください

0

$0$ 前

π / 2

$\pi/2$ ..。この間隔で線形に分布した10個のポイントを選択しましょう（選択は最適ではありません）。

x_{0} = 0

$x_0=0$ 、

x_{1} = π / 20

$x_1=\pi/20$ 、

x_{2} = 2 π / 20

$x_2=2\pi/20$ 、

x_{i} = i \cdot π / 20

$x_i=i\cdot \pi/20$ ..。各ポイントについて、このポイントでのサインとその派生物を含むプレートを計算します。これで、値を取得するときに関数を変更できるようになりました

x

$x$ 関数は最も近い値を取ります

x_{i}

$x_i$ ポイントの周りに一列に配置します

x_{i}

$x_i$ 、ゼロ付近ではない（

Δ x = x - x_{i}

$\Delta x = x - x_i$ ）。

三角測量変換の使用

さらに戻って、学校の上級クラスに戻ると、1つの非常に重要な式を思い出すことができます。

\sin (x_{0} + Δ x) = \sin (x_{0}) \cos (Δ x) + \cos (x_{0}) \sin (Δ x)

$\sin(x_0+\Delta x)=\sin(x_0)\cos(\Delta x)+\cos(x_0)\sin(\Delta x)$

そして、すべてが前の段落と同じです。間隔内のポイントを選択し、それらのサインとコサインを計算し、サイン関数を呼び出すときに、最も近いポイントを探し、上記の式を使用して、小さな値を使用してサインを計算します。

Δ x

$\Delta x$ ..。

これらの2つの方法のどちらを選択するのが良いかを考えてください。しかし、今のところ、数学から実際の計算に移ります。

浮動小数点世界における乗算の分布特性

私はインターネットにアドバイスを求めなければなりませんでした、それは何と呼ばれていますか

a \cdot (b + c) = a \cdot b + a \cdot c

$a\cdot(b+c)=a\cdot b + a\cdot c$ ..。それは配布プロパティであることが判明しました。最初のパートの終わりに私が尋ねた質問に戻りましょう。つまり、なぜ数学的に同等の式

a \cdot (b + c)

$a\cdot(b+c)$ そして

a \cdot b + a \cdot c

$a\cdot b + a\cdot c$ 浮動小数点計算で異なる結果を与えることができますか？これを説明する最も簡単な方法は、例を使用することです。4桁の精度で10進形式の浮動小数点数を処理する架空のシステムを考えてみましょう。ふりをしましょう

b = 1.000 E 0

$b=1.000\text{E}0$ 、

c = 1.234 E- 2

$c=1.234\text{E-}2$ 、および

a = 5.678 E- 2

$a=5.678\text{E-}2$ ..。まず、括弧付きの式を取り、それを段階的に計算してみましょう。各段階で四捨五入することを忘れないでください

。1）

b + c = 1.000 E 0 + 1.234 E- 2 = 1.012 E 0

$b+c=1.000\text{E}0 + 1.234\text{E-}2=1.012\text{E}0$

2）

a \cdot (b + c) = 5.678 E- 2 \times 1.012 E 0 = 5.746 E- 2

$a\cdot(b+c)=5.678\text{E-}2\times 1.012\text{E}0=5.746\text{E-}2$

受信した応答

y = 5.746 E- 2

$y=5.746\text{E-}2$

次に、同じ方法で2番目の式を段階的に計算してみましょう。

a \cdot b + a \cdot c = 5.678 E- 2 \times 1.000 E 0 + 5.678 E- 2 \times 1.234 E- 2 = 5.678 E- 2 + 7.007 E- 4 = 5.748 E- 2

$a\cdot b + a\cdot c=5.678\text{E-}2\times 1.000\text{E}0 + 5.678\text{E-}2\times 1.234\text{E-}2=5.678\text{E-}2 + 7.007\text{E-}4=5.748\text{E-}2$

受信した応答

y = 5.748 E- 2

$y=5.748\text{E-}2$

本当の答えは0.0574806652です。

ご覧のとおり、2番目のケースで得られた答えは、最初のケースよりも実際の答えにはるかに近いものです。これを指で説明すると、最初のケースで1.0に数値を追加すると想像してみてください。

c = 1.234 E- 2

$c=1.234\text{E-}2$ 最後の2桁を捨てるだけです。彼らはもうありません。2番目のケースでは、破棄は乗算後の最後に発生します。それら。2番目のケースでは、乗算演算がより正確です。

これで終わりのようですが、最初の方法を詳しく見て、計算結果がどうなるか教えてください

b + c - b

$b+c-b$ ..。そして...浮動小数点数を丸める方法があります！この例をお見逃しなく。それを理解する時間を自分に与えてください。数値の丸めは、この記事と次の記事でさらに集中的に使用されます。

この式のもう1つの特徴に注目しましょう。変数の4桁の精度では不十分だと想像してみてください。何をすべきか？そして、ここにすでに答えがあります-フォームで番号を表すために

b + c

$b+c$ 2桁の合計としてメモリに保存します。そして、それに応じて、両方の項に対して別々に操作（たとえば乗算）を実行します。この手法については、「精度を損なうことなく2つの浮動小数点数を追加する」の記事で詳しく説明しています。

前回の記事で、私もその方法を書いた

a \cdot (b + c)

$a\cdot(b+c)$ 不快な機能が1つあります。そして、それは以下の通りです。数

c

$c$ 常に数字の最後の重要な桁で切り捨てられます

b

$b$ ..。これは、数に関係なく

c

$c$ 、もし

b + c \neq b

$b+c \neq b$ 、その後、最後の記号のエラーは、小さい場合でも常に可能です

c

$c$ ..。これは、次の章のアプローチでは許可されていません。

例としてGNUライブラリを使用してどのように機能するか

どうですか？正弦の正確な計算のために、記事の冒頭で説明した2つの方法のどちらを選択しましたか？どちらの方法を選択しても、どちらも正しいです。さらに、それらは完全に同一です。私を信じて、それをチェックしてください。以下では、学校の公式を使用します。説明が簡単です。

前の記事とこの記事で得た知識を身につければ、標準ライブラリのコードを簡単に理解できます。s_sin.cファイルを開いて、そこで__sin関数を見つけましょう。

そのコードは非常に単純です。入力変数の制限に応じて異なる関数のセットを呼び出すことは容易に理解できます。この記事では、角度2 ^ -26 <| x | <0.855469の218-224コードセクションについて説明します。コードのこのセクションでは、do_sin（x、0）関数が呼び出されていることがわかります。この機能について詳しく説明します。

, dx=0 .
129-130 , abs(x)<0.126, .. x , . , , , .
136-137. , . x 2 . u x. , 0.345678. u=0.34, 0.005678.
140-142. ( s ) ( c ) x . , cos(x)=1-c, 1.0, (. ), .
143. u. , u=0.34 34. sin(u)=sn+ssn, cos(u)=cs+ccs. sn cs — «» u, ssn ccs — .
144-145. sin(u+x)=(sn+ssn)*(1-c)+(cs+ccs)*s. , , 144-145. — .

実際、私はこの方法で正弦を計算する最も簡単な部分だけを説明しました。取り残された数学はたくさんあります。たとえば、テーブルとその中の要素のサイズをどのように計算しますか？魔法の数0.126と0.855469はどこから来たのですか？テイラーの数で計算を切り落とすのはいつですか？結果を改善するためのテイラーシリーズの係数の修正。

もちろん、これはすべて興味深いものですが、客観的には、提示された方法には多くの欠点があります。サイン（s）とコサイン（c）を同時に計算する必要があり、テイラーシリーズ^1の2倍の計算が必要です。ご覧のとおり、表形式の値による乗算も無料ではありません。また、3520バイトのテーブルをRAMに格納することはもちろん問題ではありませんが、（キャッシュ内であっても）それにアクセスするにはコストがかかる可能性があります。

したがって、次のパートでは、プレートを取り除き、間隔[0.126、0.855469]の正弦を直接計算しますが、最初の章よりも正確に計算します。

終了する前に-クイックウィットの質問。この例では、数の大きなは52776558133248 = 3 * 2である⁴⁴。そのような数は、たとえば2 ⁴⁵ではなく、どこから来たのですか？質問をより正確に定式化します。たとえば、2 ^{N + 1}ではなく、数値を丸めるときに、数値3 * 2 ^Nが最適なのはなぜですか？別の質問ですが、数値を整数に丸めるためにどのNを選択する必要がありますか？¹サインとコサインを同じ角度から同時に計算すると、このアプローチの大きな利点が現れる可能性があることに注意してください。2番目の関数はほとんど無料で計算できます。

正弦関数の例を使用した浮動小数点数の正確で高速な計算。パート2：libm

もう一度いくつかの数学

与えられた値の周りのシリーズ展開

三角測量変換の使用

浮動小数点世界における乗算の分布特性

例としてGNUライブラリを使用してどのように機能するか

More articles: