現代のコンパイラで使用されている数学ライブラリが正しく丸められた結果をもたらさないことがあることを気にする人がほとんどいないかのように、ロシア語でこの問題に関する情報を見つけるのが難しいことに驚きました。私はそのような数学的ライブラリの開発に取り組んでいるだけなので、この状況が心配です。外国の文献では、この問題は十分にカバーされているので、西洋の情報源とまだ少し個人的な経験に頼って、人気のある科学的な形でロシア語でそれを提示することにしました。

---

友だち、あなたの便宜のために、記事はビデオプレゼンテーション形式（約34分）でも利用できます。この形式は、プレゼンテーションに多くの説明資料があるため、頭の中で必要な数学的画像を作成するのが難しい読者に適しています。ビデオの情報は、記事の内容と完全に同じです。ご都合の良いときに行動してください。

---

これは科学的なものではなく、人気のある科学的な記事であることを繰り返します。これを読んだ後、簡単にこれを知ることができます。

• 浮動小数点演算で機能する超越的な基本関数（exp、sin、log、coshなど）は正しく丸められず、最後のビットでエラーが発生することがあります。
• エラーの理由は、必ずしも開発者の怠惰や資格の低さにあるとは限りませんが、現代科学がまだ克服できていない1つの基本的な状況にあります。
• «», - .
• , , , , exp2(x) pow(2.0, x).

この記事を理解するには、IEEE-754フローティングポイント形式に精通している必要があります。たとえば、これが次のとおりであることを少なくとも理解していれば十分です。0x400921FB54442D18-2倍精度形式（binary64、またはdouble）の数値pi、つまり、このレコードの意味を理解しているだけです。私はそのような変換をその場で行うことができることを要求しません。そして、この記事の丸めモードについて思い出させます。これはストーリーの重要な部分です。西洋の文学からの用語や引用があるので、「プログラマー」の英語を知っていることも望ましいですが、あなたはオンライン翻訳者と一緒にやっていくことができます。



最初に例を挙げて、会話の主題が何であるかをすぐに理解できるようにします。ここでコードをC ++で示しますが、これがあなたの言語でない場合でも、書かれていることを簡単に理解できると確信しています。このコードを見てください：

#include <stdio.h>
#include <cmath>

int main() {
  float x = 0.00296957581304013729095458984375f;  // ,  .
  float z;
  z = exp2f(x);  // z = 2**x  .
  printf ("%.8f\n", z);  //      8   .
  z = powf(2.0f, x);  // z = 2**x  
  printf ("%.8f\n", z);  //   .
  return 0;
}

数値xは、floatタイプで正確に表現できるように、つまり、コンパイラーが丸めずにバイナリコードに変換するように、意図的にこのような有効桁数で記述されています。結局のところ、一部のコンパイラはエラーなしで丸めることができないことをよく知っています（わからない場合は、コメントで示してください。例を含む別の記事を作成します）。次にプログラムで^2xを計算する必要がありますが、2つの方法で計算しましょう。関数exp2f（x）と、2つのpowf（2.0f、x）の明示的な指数です。もちろん、結果は異なります。これは、基本機能がすべての場合に正しく機能するとは限らないことを前述したため、これを示すために特別に例を選択しました。出力は次のとおりです。

1.00206053
1.00206041

Microsoft C ++（19.00.23026）、Intel C ++ 15.0、GCC（6.3.0）、Clang（3.7.0）の4つのコンパイラーからこれらの値が得られました。それらは1つの最下位ビットが異なります。これらの数値の16進コードは次のとおりです。

0x3F804385  // 
0x3F804384  // 

この例を覚えておいてください。少し後で問題の本質を見ていきますが、今のところ、より明確な印象を得るために、他のいくつかの基本関数を使用した倍精度データタイプ（double、binary64）の例を参照してください。結果を表に示します。正解（利用可能な場合）の最後に*が付いています。







私がわざと見つけられないような完全にユニークなテストを行ったという印象を受けないでください。もしそうなら、floatデータタイプの^2x関数のすべての可能な分数引数の完全な列挙をひざまずいてみましょう。他の引数は指数フィールドの値のみが異なり、関心のない結果をもたらすため、0から1の間のxの値のみに関心があることは明らかです。あなた自身が理解しています：

$$

$$2^x = 2^{[x]}\cdot2^{\{x\}}.$$

そのようなプログラムを書いた後（非表示のテキストは以下にあります）、exp2f関数と、0から1までの間隔xで生成される誤った値の数を確認しました。







以下のプログラムから、テストされた引数xの数が197612997であったことが明らかです。たとえば、Microsoft C ++は、それらのほぼ1パーセントについて^2x関数を誤って計算していることがわかりました。GCCとClangのファンの皆さん、喜ばないでください。この関数がこれらのコンパイラで正しく実装されているだけですが、他のコンパイラではエラーがいっぱいです。





私はこれらの例で読者を退屈させません。ここでの主なことは、超越関数の最新の実装が最後のビットを誤って丸めることができ、さまざまなコンパイラがさまざまな場所で間違いを犯す可能性があることを示すことでしたが、どれも正しく機能しません。ちなみに、IEEE-754標準では最後のビットでこのエラーが許可されていますが（これについては後で説明します）、それでも私には奇妙に思えます。これは大きなデータタイプですが、フロートはブルートフォースでチェックできます。そんなに大変でしたか？まったく難しいことではありません、そして私はすでに例を示しました。



私たちの列挙コードには、正しい計算の「自己記述」関数が含まれています2 ^x10次​​の近似多項式を使用し、そのような多項式は、たとえばMapleコンピュータ代数システムで自動的に導出されるため、数分で記述されました。多項式の条件を設定して54ビットの精度を提供するだけで十分です（この関数の場合、2 ^x）。なぜ54？しかし、問題の本質を説明し、理論的にはこの問題を攻撃する試みはすでに行われていますが、原則として、4倍精度のデータタイプ（binary128）に対して高速で正確な超越関数を作成することが不可能である理由を説明した直後にすぐにわかります。

デフォルトの丸めとその問題

数学ライブラリの開発に没頭していない場合は、IEEE-754標準に準拠した浮動小数点数のデフォルトの丸め規則を忘れても問題はありません。したがって、私はあなたにそれを思い出させます。すべてをよく覚えているなら、とにかく少なくともこのセクションの終わりを見てください、あなたは驚きに満ちています：私はあなたに数を切り上げることが非常に難しいかもしれない状況をあなたに見せます。



「切り上げ」（プラスの無限大へ）、「切り下げ」（マイナスの無限大へ）、「ゼロへの丸め」を名前で簡単に思い出すことができます（どちらかといえば、ウィキペディアがあります））。プログラマーにとっての主な問題は、「最も近いものに丸めるが、最も近いものから等しい距離の場合、最後の桁が偶数であるものに」丸めることで発生します。はい、これはこの丸めモードがどのように変換されるかであり、西洋の文献は要するに「最も近いものを偶数に丸める」と呼んでいます。



この丸めモードはデフォルトで使用され、次のように機能します。計算の結果、マンティッサの長さが結果のデータタイプに対応できる長さよりも長いことが判明した場合、2つの可能な値の最も近い値に丸めが実行されます。ただし、元の数値が最も近い2つの数値のちょうど中間にあることが判明した場合、最後のビット（丸め後）が偶数、つまりゼロに等しい結果が選択される場合があります。バイナリ小数点の後に2ビットに丸める必要がある4つの例を考えてみましょう。

1. 1.00 1 001を丸めます。小数点以下の3番目のビットは1ですが、別の6番目のビットである1があります。これは、元の数値が1.00よりも1.01に近いため、丸めが上がることを意味します。
2. 1,001000. , 1,00 1,01, .
3. 1,011000. 1,01 1,10. , .
4. 1,010111. , 1,01, 1,10.

これらの例から、すべてが単純に見えるかもしれませんが、そうではありません。実際には、2つの値の中間にあるかどうかを確実に判断できない場合があります。例を参照してください。我々は再び小数点以下2ビットにラウンドしたいとし：



1.00 1 0000000000000000000000000000000000001



丸めが数1.01に、あること、アップであることを今、あなたには自明です。ただし、小数点以下40ビットの数値を見ています。アルゴリズムが40ビットの精度を提供できず、30ビットしか達成できない場合はどうなりますか？それからそれは別の数を与えるでしょう：



1.00 1 000000000000000000000000000



40番目の位置（アルゴリズムでは計算できません）に大切なものがあることに気づかずに、この数値を切り捨てて1.00を取得しますが、これは間違っています。あなたは最後のビットを間違って切り上げました-それが私たちの議論の主題です。上記のことから、2番目のビットだけを正しくするためには、最大40ビットの関数を計算する必要があることがわかります。うわー！そして、ゼロの「機関車」がさらに長くなることが判明した場合はどうなりますか？これについては、次のセクションで説明します。



ちなみに、これは多くのコンパイラが浮動小数点数の10進表記を結果のバイナリ形式に変換するときに犯す間違いです。プログラムコードの元の10進数が、正確に表現できる2つのバイナリ値の中間に近すぎると、正しく丸められません。しかし、これはこの記事のトピックではなく、別の話の理由です。

最後の重要なビットを丸める問題の本質

この問題は2つの理由で現れます。1つ目は、時間のかかる計算を意図的に拒否し、速度を優先することです。この場合、指定された精度が守られている限り、応答にどのビットが含まれるかは二次的な問題です。2番目の理由は、私たちの会話の主要な主題であるテーブルメーカーのジレンマです。両方の理由をさらに詳しく考えてみましょう。

最初の理由

もちろん、超越関数の計算は、いくつかの近似方法、たとえば、多項式を近似する方法、または（まれに）系列展開によって実装されることを理解しています。計算をできるだけ早く行うために、開発者は、アルゴリズムがマンティッサの最後のビットの値の半分を超えないエラーを許可する限り、数値メソッドの反復をできるだけ少なく実行する（または可能な限り最小の多項式を取る）ことに同意します。文献では、これは0.5ulpと書かれています（ulp =最後の単位）。



たとえば、間隔（0.5; 1）内のfloat型の数xについて話している場合、値ulp = ^2-23です。間隔（1; 2）でulp = ^2-22。言い換えると、xが間隔（0; 1）にある場合、2 ^xは間隔（1,2）にあり、0.5ulpの精度を確保するには、大まかに言えば、EPS = 2 ^-23を選択する必要があります（したがって、「エラー」または「精度」と呼ばれる一般の人々では、定数「イプシロン」を示します。あなたが好きなように、欠点を見つけないでください）。



適用される計算の場合、これで十分ですが、最後のビットが絶対結果と一致しない可能性があるという事実は、ほとんどのプログラマーにとって重要ではありません。ビットが何であるかではなく、精度がどうなるかが重要だからです。



わからない方のために、10進数制の例をあげます。 1.999999と2.0の2つの数値があります。 1つ目はプログラマーが受け取ったものであり、2つ目は無限の可能性があった場合に得られるべきものの標準であるとしましょう。それらの違いはわずか100万分の1です。つまり、回答はEPS = ^10-6で計算されました。ただし、この回答には正しい番号は1つではありません。悪いですか？いいえ、アプリケーションプログラムの観点からは、これは紫色です。プログラマーは答えを小数点以下2桁に切り上げて、2.00（たとえば、通貨についての$ 2.00）を取得します。これ以上は必要ありませんが、彼が私のプログラムにEPS = 10 ^-6を入れて、それでうまくいき、中間計算のエラーのマージンを取り、問題を正しく解決しました。



言い換えれば、混同しないでください。正しいビット（または桁）の精度と数は2つの異なるものです。正確さが必要な人（これはほぼ100％のプログラマーです）、議論された問題はそれらにはまったく関係ありません。正しく丸められた参照に一致するためにビットシーケンスを必要とする人は誰でも、たとえば、基本関数のライブラリの開発者など、この問題について非常に心配しています。それにもかかわらず、一般的な開発のためにこれについて知っていることは誰にとっても有用です。



これが問題の最初の方向であったことを思い出させてください。これは意図的な解決策であるため、答えの最後の部分が間違っている可能性があります。主なことは、0.5ulp（またはそれ以上）の精度を維持することです。したがって、数値アルゴリズムは、それが非常に高速に機能する場合にのみ、この条件からのみ選択されます。同時に、この規格では、最後のビットを正しく丸めることなく、基本関数を実装できます。私は[1、セクション12.1]（英語）を引用します：

浮動小数点演算に関するIEEE754標準の1985バージョンでは、基本機能に関して何も指定されていませんでした。これは、正しく丸められた関数は、少なくとも一部の入力引数に対しては遅すぎると長年信じられてきたためです。それ以来状況は変化し、2008年版の規格では、一部の関数を正しく丸めることを推奨しています（ただし必須ではありません）。

以下は、推奨されているが正しく丸める必要のない関数です。

2番目の理由

最後に、会話のトピックに到達しました：テーブルメーカーのジレンマ（TMDと略されます）。その名前をロシア語に適切に翻訳することはできませんでした。それはWilliamKahan（IEEE-754の創設者）によって記事[2]で紹介されました。おそらく記事を読めば、その名前がまさにそれである理由を理解するでしょう。要するに、ジレンマの本質は、完全に計算された結果zの無限ビットレコードを自由に使用できるかのように、関数z = f（x）を完全に正確に丸める必要があるということです。しかし、無限のシーケンスを取得できないことは誰にとっても明らかです。それでは何ビットかかりますか？上記では、丸め後に少なくとも2つの正しいビットを取得するために、結果の40ビットを確認する必要がある場合の例を示しました。そして、TMD問題の本質は、事前にわからないということです。、必要な数の丸め後に正しいビット数を取得するために、zの値を計算するビット数まで。百または千がある場合はどうなりますか？事前にはわかりません！



たとえば、前述したように、関数2 ^xの場合、マンティッサの分数部分が23ビットしかないデータ型floatの場合、例外なくすべての可能なx引数に対して正しく丸めが行われるように、^2-54の精度で計算を実行する必要があります。徹底的な検索でこの推定値を取得することは難しくありませんが、他のほとんどの関数、特にタイプdoubleまたはlong double（それが何であるかを知っている場合は「クラス」と入力）の場合、そのような推定値は不明です。



なぜこれが起こっているのかをすでに理解しましょう。この記事の最初の例としてfloatデータタイプを意図的に示し、覚えておいてもらいました。このタイプでは32ビットしかないため、見やすくなります。他のデータタイプでも状況は似ています。



数値x = 0.00296957581304013729095458984375から始めました。これは、floatデータタイプで正確に表現できる数値です。つまり、丸めることなくバイナリfloatシステムに変換できるように記述されています。 2 ^xを計算し、無限の精度の計算機がある場合は、次のようになります（確認できるように、計算はオンラインのWolframAlphaシステムで行われます）：



1.0020604729652405753669743044108123031635398201893943954577320057 .. ..



：レッツの64ビットは十分でしょう言わせて、バイナリにこの番号を変換



1.00000000100001110000100 1 000000000000000000000000000001101111101



丸めビット（小数点以下24ビット）は下線が引かれています。質問：どこを丸めるのですか？上か下？明らかに、あなたは十分なビットを見て決定を下すことができるので、これを知っています。しかし、注意深く見てください...



丸めビットの後、29個のゼロがあります。これは、最も近い2つの数値の中間に非常に近く、丸めの方向が変わるため、少し下に移動するだけで十分であることを意味します。しかし、問題は、このシフトはどこにあるのかということです。数値アルゴリズムは、段階的に、さまざまな側面から正確な値に近づくことができます。これらの29個のゼロをすべて通過し、この「機関車」の最後のゼロの値を超える精度に達するまで、丸めの方向はわかりません。 ..。実際、正解が次のようになっている場合はどうなりますか：



1.00000000100001110000100 0 11111111111111111111111111111？



その後、丸めがダウンします。



精度が小数点以下54ビットに達するまでこれはわかりません。 54番目のビットが正確にわかっている場合、2つの最も近い数値のどちらが実際に近いかを正確に知ることができます。そのような番号が呼び出される最も困難なツー丸点[1、セクション12.3（丸めのための重要な点）、及び数54が呼び出され、硬度対ラウンド、および引用した書籍の文字Mで示されています。



丸めの複雑さ（m）は、特定の関数f（x）のすべての引数について、および事前に選択された範囲について、関数f（x）が最後のビットに正しく丸められるようにするために必要な最小ビット数です（丸めモードが異なると、異なる場合があります）値m）。言い換えると、データ型floatの場合、および丸めモードの範囲（0; 1）から「最も近い偶数」の丸め時間m = 54までの引数xの場合。これは、間隔（0; 1）からの絶対にすべてのxについて、同じ精度のESP = 2 ^-54をアルゴリズムに入れることができ、すべての結果が2進小数点以下23ビットに正しく丸められることを意味します。



実際、一部のアルゴリズムは正確な結果を提供でき、53ビットおよび52ビットに基づいて、ブルートフォースがこれを示しますが、理論的には正確に54が必要です。ブルートフォースをクランクアウトする可能性がなければ、「チート」することはできません。上記のブルートフォースプログラムで行ったように、2、3ビット節約します。必要以上に低い次数の多項式を取りましたが、運が良かったという理由だけで、それでも機能します。



したがって、丸めモードに関係なく、2つの状況が考えられます。丸め領域でゼロの「蒸気機関車」が発生するか、1の「蒸気機関車」が発生します。超越関数f（x）を計算するための正しいアルゴリズムのタスクは、精度がこの「蒸気機関車」の最後のビットの値を超えるまで、そしてfを計算するための数値アルゴリズムのその後の変動の結果として正確に明らかになるまで、この関数の値を改良することです。 （x）ゼロは1になりません、またはその逆です。すべてが安定し、アルゴリズムが「蒸気機関車」の限界を超える精度に達するとすぐに、ビット数が無限であるかのように丸めることができます。そして、この丸めは正しい最後のビットで行われます。しかし、これはどのように達成できますか？

「クラッチ」

前述のように、主な問題は、丸めビットの直後に来る0または1の機関車を克服するためのアルゴリズムを取得することです。機関車が克服され、全体として見ると、これは、これらの0または1がすでに正確に計算されており、丸めがどの方向で発生するかを正確に知っているという事実と同等です。しかし、機関車の長さがわからない場合、どのようにアルゴリズムを設計できますか？

最初の「クラッチ」

答えは明白であるように読者には思われるかもしれません。無限の精度で算術を取り、意図的に過剰なビット数を入れ、それが十分でない場合は、別のビットを入れて再計算します。一般的に、それは正しいです。これは、コンピューターの速度とリソースが特別な役割を果たさない場合に行われます。このアプローチには名前があります：Zivのマルチレベル戦略[1、セクション12.3]。その本質は非常に単純です。アルゴリズムは、いくつかのレベルでの計算をサポートする必要があります。迅速な予備計算（ほとんどの場合、最終的なものであることが判明）、低速ですがより正確な計算（最も重要な場合に保存）、さらに低速ですが、さらに正確な計算（絶対に「悪い場合」） 「しなければならなかった）など。



圧倒的多数の場合、0.5ulpより少し高い精度で十分ですが、「機関車」が出てきたら上げます。 「蒸気機関車」が残っている限り、数値法のさらなる変動がこの「蒸気機関車」に影響を与えないことが厳密に明らかになるまで、精度を上げます。したがって、たとえば、この場合、ESP = 2 ^-54に達した場合、54番目の位置にユニットが表示されます。これは、いわば、機関車をゼロから「保護」し、^2-53以上の値の減算が発生しないことを保証します。ゼロは1にならず、丸めビットをゼロにドラッグします。



これは人気のある科学プレゼンテーションであり、Zivの丸めテストと同じです。このテストでは、[1、第12章]、または[3、セクション]で、必要な精度を達成したかどうかを1つのステップで確認する速度が示されています。 10.5]。



このアプローチの問題は明らかです。各超越関数f（x）を計算するためのアルゴリズムを設計する必要があります。これにより、ピースの過程で計算の精度を高めることができます。ソフトウェアの実装の場合、これはまだそれほど怖いものではありません。たとえば、ニュートンの方法では、大まかに言えば、各反復で小数点以下の正確なビット数を2倍にすることができます。 「十分」になるまで2倍にすることができますが、これはかなり時間のかかるプロセスですが、逆関数f ^-1を計算する必要があるため、ニュートンの方法が常に正当化されるとは限らないことを認めなければなりません。（x）、場合によっては、f（x）自体を計算するよりも簡単ではない場合があります。ハードウェアの実装では、「Ziva戦略」は完全に不適切です。プロセッサに組み込まれたアルゴリズムは、すでに事前設定されたビット数で一連のアクションを実行する必要があります。この数が事前にわからない場合、これを実装するのは非常に問題があります。株式を取ります？いくら？



問題を解決するための確率論的アプローチ[1、セクション12.6]により、mの値を推定することができます（これはビット数であり、正しい丸めに十分であることを忘れないでください）。確率的な意味での「機関車」の長さは、数のマンティッサの長さよりわずかに長いことがわかります。したがって、ほとんどの場合、mはマンティッサの値の2倍強を取るだけで十分であり、ごくまれなケースでのみ、さらに多く取る必要があります。私はこの作品の作者を引用します：「実際には、mは2pよりわずかに大きくなければならない」（彼らはpを持っています-整数部分と一緒にマンティッサの長さ、つまりフロートの場合はp = 24）。さらに本文では、そのような戦略でのエラーの確率はゼロに近いが、それでも正であることが示され、これは実験によって確認されています。



それでもなお、mの値をさらにとらなければならない場合があり、最悪の場合は事前にわかっていません。最悪の場合の理論的な推定値は存在しますが[1、セクション12.7.2]、それらは考えられない数百万のビットを生成します。これは良くありません。これは引用された作業からの表です（これは-ln（2）からln（2）までの間隔での関数exp（x）用です）：

2番目の「クラッチ」

実際には、mはそれほど大きくはありません。そして、最悪のケースを決定するために、「徹底的な事前計算」と呼ばれる2番目の「クラッチ」が適用されます。データタイプfloat（32ビット）の場合、関数fに1つの引数（x）があると、xのすべての可能な値を簡単に「実行」できます。この問題は、複数の引数（pow（x、y）の中で）を持つ関数でのみ発生します。そのため、そのようなことは考えられませんでした。 xのすべての可能な値を確認した後、各関数f（x）および各丸めモードの定数mを計算します。次に、ハードウェアに実装する必要のある計算アルゴリズムは、2- ^mの精度を提供するように設計されています。そうすれば、f（x）の丸めはすべての場合に正しいことが保証されます。



ダブルタイプ（64ビット）の場合、単純な列挙はほとんど不可能です。しかし、彼らは整理しています！しかし、どのように？答えは[4]にあります。それについて簡単に説明します。



関数f（x）のドメインは非常に小さなセグメントに分割されているため、各セグメント内でf（x）をb-ax形式の線形関数に置き換えることができます（もちろん、係数aとbはセグメントごとに異なります）。これらのセグメントのサイズは分析的に計算されるため、このような線形関数は実際には各セグメントの元の関数とほとんど区別できません。



次に、スケーリングとシフトの操作を行った後、次の問題が発生します。直線のb-axを整数点に「十分に近づける」ことができるでしょうか。



はいまたはいいえの答えを出すのは比較的簡単であることがわかります。つまり、潜在的に危険なポイントが直線に近い場合は「はい」、原則としてそのようなポイントが直線に近づかない場合は「いいえ」です。この方法の利点は、実際には圧倒的多数のケースで「いいえ」の回答が得られ、めったに得られない「はい」の回答が、どの特定のポイントが重要であるかを判断するために徹底的な検索でセグメントを通過することを強制することです。



それでも、f（x）の引数の繰り返しは何度も削減され、double（binary64）やlong double（80ビット！）などの数値のブレークポイントを検出できます。これは、スーパーコンピューター、そしてもちろんビデオカードで行われます...マイニングからの自由な時間に。ただし、binary128データタイプをどうするかはまだ誰にもわかりません。そのような数のマンティッサの分数部分は112ビットであることを思い出させてください。したがって、この主題に関する外国の文献では、「私たちは希望する...」（「私たちは希望する...」）で始まる半哲学的な推論しか見つけることができません。



整数点の近くの線の通過をすばやく決定できる方法の詳細は、ここでは不適切です。プロセスを学びたい方のために、私は直線とZの間の距離見つける問題に向けて、より慎重に探してお勧め^2を記事に、例えば、[5]。改良されたアルゴリズムについて説明します。これは、構築の過程で、最大の共通除数を見つけるための有名なユークリッドのアルゴリズムに似ています。 [4]と[5]から同じ図を示します。







これは、問題のさらなる変換を示しています。超越関数ごとに異なる間隔で丸める最悪のケースを含む広範なテーブルがあります。それらは[1セクション12.8.4]と[3、セクション10.5.3.2]、および別の記事、たとえば[6]にあります。



そのようなテーブルからランダムな行を取得して、いくつかの例を示します。これらはすべてのxにとって最悪のケースではないことを強調しますが、ごく一部の間隔でのみです。興味がある場合はソースを参照してください。







表の読み方は？値xは、16進浮動小数点の二重表記で指定されます。最初に、予想どおり、先頭に1つあり、次にマンティッサの小数部分の52ビットと文字Pがあります。この文字は、「2の累乗」の後に度が続くことを意味します。たとえば、P-23は、指定された^{マンティッサに2-23を}掛ける必要があることを意味します。



さらに、関数f（x）が無限の精度で計算され、最初の53ビットが（丸めなしで）切り取られると想像してください。 f（x）列に示されているのは、これらの53ビット（そのうちの1つはコンマまで）です。後続のビットは最後の列に示されています。最後の列のビットシーケンスの「次数」記号は、ビットの繰り返し回数を意味します。つまり、たとえば10 ⁵³1011は、最初に1に等しいビット、次に53個のゼロ、次に1011になることを意味します。次に、省略形。これは、一般に、残りのビットをまったく必要としないことを意味します。



さらに、それは技術の問題です-別々に取られた関数の各間隔の最悪のケースを知っており、この間隔に対して、最悪のケースを正確にカバーするような近似を選択できます。わずか数年のスーパーコンピューターコンピューティングで、基本機能の高速で正確なハードウェア実装を作成することが可能です。問題は小さいです。少なくともコンパイラ開発者にこれらのテーブルを使用するように教えることは残っています。

なぜこれが必要なのですか？

素晴らしい質問です！結局のところ、プログラマーのほぼ100％は、正しく丸められた最後のビット内で基本関数を知る必要がない（多くの場合、ビットの半分も必要ない）と繰り返し話しましたが、なぜ科学者はスーパーコンピューターを駆動し、テーブルをコンパイルして「役に立たない」問題を解決するのですか？



まず、課題は基本です。正確な丸めのために正確な丸めを取得しないことはかなり興味深いですが、原則として、この興味深い問題をどのように解決できるかを理解するために、その解決策が私たちに明らかにする計算数学の秘密は何ですか？これらの秘密を他のタスクでどのように使用できますか？基礎科学-そのようなものです。何十年もの間、ある種の「ナンセンス」を行うことができます。その後、100年後、この「ナンセンス」のおかげで、他の分野で科学的な進歩が起こります。



第二に、コードの移植性の問題。関数が結果の最後のビットを必要な方法で処理できる場合、それは、異なるプラットフォームおよび異なるコンパイラーで、わずかに異なる結果が得られる可能性があることを意味します（指定されたエラー内であっても）。これは重要でない場合もありますが、特にプログラムに1つのプラットフォームに表示されるエラーが含まれているが、結果のビットが異なるために別のプラットフォームには表示されない場合は、重要な場合があります。しかし、なぜ私はあなたに異なるプログラムの振る舞いに関連するよく知られた頭痛を説明しているのですか？あなたは私なしでこれをすべて知っています。どのようにコンパイルされていても、すべてのプラットフォームでまったく同じように機能する数学システムがあれば素晴らしいと思います。このためにあなたは正しくする必要があります 最後のビットを丸めます。

ソースのリスト

[1] Jean-Michel Muller、「基本関数：アルゴリズムと実装」、2016年



[2] William Kahan、「対数が半分賢すぎる」、2004年



[3] Jean-Michel Muller、「浮動小数点演算のハンドブック」 、2018



[4]VincentLefèvre、Jean-Michel Muller、「正しく丸められた超越に向けて」、IEEE TRANSACTIONS ON COMPUTERS、VOL。47、いいえ。1998年11月11日。pp。1235-1243



[5] VincentLefèvre 。「セグメントとZとの間の距離に新しい結果²」。正確な丸めへの適用。コンピュータ算術に関する第17回IEEEシンポジウム-Arith'17、2005年6月、ケープコッド、マサチューセッツ州、

アメリカ合衆国。pp.68-75



[6]VincentLefèvre、Jean-Michel Muller、「倍精度での基本関数の正しい丸めの最悪のケース」、Rapport de recherche（INSTITUT NA TIONAL DE RECHERCHE EN INFORMA TIQUE ET EN AUTOMA TIQUE）n˚4044-2000年11月- 19ページ。