日常の算術演算に関しては、計算の有限精度の問題はそれほど恐ろしく見えません。 そして、結果が正しく得られることを確認する最善の方法は、異なる精度で得られた値を比較することです。
たとえば、単精度と倍精度で得られた計算が一致する場合、少なくとも単精度に匹敵する精度で、結果に自信が生まれます。 ここでは、 比較的単純な算術問題であっても、 数値の表現における可変精度のこのような安定性が、そのような信頼性の基礎として役立たないことを実証する興味深い例を挙げたいと思います。
引数の特定の値（以下に示す）に対して2つの変数の関数の値を計算する必要があります。
この例は、 C-XSCライブラリ（（主に）任意の精度の間隔科学計算のためのC言語クラスシステム）を扱っていたときに気づきました。 このライブラリは、さまざまな数値アルゴリズムの計算安定性に関する実用的な研究に最適です。
Pythonで浮動小数点計算をエミュレートするには、 mpmathパッケージをインストールします。 原則として、IEEE 754の精度に制限する場合、NumPyを使用できますが、IEEE 754のフレームワークで得られた結果が正しくないことを示す必要があります。
関数f （ a 、 b ）の値は、 a = 77617およびb = 33096で計算します。
# coding: utf-8 from mpmath import * def f(a,b): ''' From: Ramp SM Algorithms for verified inclusions -- theory and practice, USA, NY, 1988. ''' return (333.75 - a * a) * b ** 6 + a * a * (11 * a * a * b * b - 121 * b ** 4 - 2) + 5.5 * b ** 8 + a/(2.0*b) # mp.dps = 8 a = mpf(77617) b = mpf(33096) print ' f(a, b) :', f(a, b) # mp.dps = 16 a = mpf(77617) b = mpf(33096) print ' f(a, b) :', f(a, b)
読者は、 mpmathでdpsを使用して精度を設定することは、実際の倍精度をエミュレートする際に適切なアプローチではないことに気付くかもしれません。 IEEEのフレームワーク内で倍精度および/または単精度について話している場合は、おそらく、バイナリシステムでそれらの特性を記述する方が適切です。 ただし、ここではそれほど重要ではありません。 しかし、 mp.dpsの使用はより簡単に解釈されます-つまり 数値の表現における10進数の有効数字の数として。
コードを実行すると、次のものが得られます。
f(a, b) : 1.1726039 f(a, b) : 1.172603940053179
とても説得力がありますよね？ 値だけが間違っています！ 与えられたaとbの正しい値は、通常1未満です！
次の計算で例を補完します。
for i in range(8, 40): mp.dps = i a = mpf(77617) b = mpf(33096) print ' dps={0}, f(a,b)={1}'.format(mp.dps, f(a,b))
取得（一部の行は省略）：
dps=8, f(a,b)=1.1726039 dps=9, f(a,b)=1.17260394 dps=10, f(a,b)=-7.737125246e+25 ... dps=13, f(a,b)=1.172603940053 dps=14, f(a,b)=1.1726039400532 dps=15, f(a,b)=1.17260394005318 dps=16, f(a,b)=1.172603940053179 dps=17, f(a,b)=-9.2233720368547758e+18 ... dps=28, f(a,b)=1.172603940053178631858834905 dps=29, f(a,b)=1.1726039400531786318588349045 dps=30, f(a,b)=1.17260394005317863185883490452 ... dps=36, f(a,b)=-0.827396059946821368141165095479816292 dps=37, f(a,b)=-0.827396059946821368141165095479816292 dps=38, f(a,b)=-0.827396059946821368141165095479816292 dps=39, f(a,b)=-0.827396059946821368141165095479816291999
安定性にもかかわらず、一部の値は通常の値とは大きく異なることがわかります。
精度dps = 36以上で得られた値は正しいとすぐに言います。 しかし、1.17260の値で見られたように、さまざまな精度での結果の安定性でさえ、その正確性を保証できないため、精度がさらに向上してもジャンプが発生しないことをどのようにして知ることができますか？
幸いなことに、この質問はmpmathパッケージを使用して回答することもできます。 このパッケージを使用すると、間隔の計算を実行できます。これにより、引数を異なる精度で提示するときに、関数の可能な値の範囲を評価できます。
mpmath区間演算装置を使用して計算を実行します。
for j in range(6, 40): iv.dps = j a = iv.mpf(77617) b = iv.mpf(33096) print '={0}: f(a, b)={1}'.format(mp.dps, f(a,b))
次のものが得られます。
dps=6: f(a, b)=[-5.0706024009e+30, 6.3382542101e+30] dps=7: f(a, b)=[-3.16912650057e+29, 3.16912654779e+29] dps=8: f(a, b)=[-5.94211218857e+28, 2.971056097974e+28] dps=9: f(a, b)=[-3.7138201178561e+27, 4.9517601582944e+27] dps=10: f(a, b)=[-1.54742504910673e+26, 3.09485009825849e+26] dps=11: f(a, b)=[-1.934281311383407e+25, 3.86856262277385e+25] dps=12: f(a, b)=[-4.8357032784585167e+24, 3.6267774588444373e+24] dps=13: f(a, b)=[-3.02231454903657294e+23, 2.26673591177745118e+23] dps=14: f(a, b)=[-2.833419889721787128e+22, 2.833419889721790484e+22] dps=15: f(a, b)=[-3.5417748621522339103e+21, 3.5417748621522344346e+21] dps=16: f(a, b)=[-442721857769029238784.0, 442721857769029246976.0] dps=17: f(a, b)=[-27670116110564327424.0, 27670116110564327456.0] dps=18: f(a, b)=[-3458764513820540927.0, 2305843009213693953.5] dps=19: f(a, b)=[-432345564227567614.828125, 288230376151711745.179688] dps=20: f(a, b)=[-36028797018963966.8274231, 18014398509481985.17260742] dps=21: f(a, b)=[-3377699720527870.8273963928, 1125899906842625.172604084] dps=22: f(a, b)=[-140737488355326.827396061271, 422212465065985.172603942454] dps=23: f(a, b)=[-17592186044414.8273960599472, 17592186044417.17260394006735] dps=24: f(a, b)=[-1099511627774.8273960599468637, 3298534883329.17260394005325] dps=25: f(a, b)=[-137438953470.827396059946822859, 412316860417.172603940053178917] dps=26: f(a, b)=[-17179869182.82739605994682137446, 17179869185.17260394005317863941] dps=27: f(a, b)=[-2147483646.8273960599468213681761, 2147483649.1726039400531786320407] dps=28: f(a, b)=[-268435454.827396059946821368142245, 268435457.172603940053178631864532] dps=29: f(a, b)=[-8388606.827396059946821368141165871, 8388609.172603940053178631858840746] dps=30: f(a, b)=[-1048574.8273960599468213681411651475, 1048577.1726039400531786318588349559] dps=31: f(a, b)=[-131070.827396059946821368141165095792, 131073.172603940053178631858834907439] dps=32: f(a, b)=[-8190.827396059946821368141165095483, 1.172603940053178631858834904520184761] dps=33: f(a, b)=[-1022.8273960599468213681411650954798445, 1.1726039400531786318588349045201837979] dps=34: f(a, b)=[-126.82739605994682136814116509547981678, 1.17260394005317863185883490452018372567] dps=35: f(a, b)=[-6.827396059946821368141165095479816292382, 1.172603940053178631858834904520183709123] dps=36: f(a, b)=[-0.82739605994682136814116509547981629200549, -0.82739605994682136814116509547981629181741] dps=37: f(a, b)=[-0.827396059946821368141165095479816292005489, -0.827396059946821368141165095479816291981979] dps=38: f(a, b)=[-0.8273960599468213681411650954798162919996113, -0.827396059946821368141165095479816291998142] dps=39: f(a, b)=[-0.82739605994682136814116509547981629199906033, -0.82739605994682136814116509547981629199887666]
関数の可能な値の間隔の変化を追跡することは興味深いです。36以上の有効な10進数の精度を使用する場合にのみ安定しますが、徐々に狭くなります。
間隔の計算から、 dps = 36以上の小数桁から正確に得られた結果のみを信頼する必要があることが非常に明確になります。
この例は、浮動小数点計算の実行がいかに慎重であるかを明確に示しており、128ビット（たとえば、PythonとNumPyについて話す場合はnumpy.float128 ）の精度でさえ不十分な場合があります。 また、さまざまな精度で得られた結果の安定性を信頼できないことも示しています。 間隔計算装置の使用は、この場合の解決策の1つであり、適切な結果を得るために必要な精度を評価することができます。