57

楽しみのために、そして本当に簡単だったので、 Grafting numbersを生成する短いプログラムを作成しましたが、浮動小数点の精度の問題により、より大きな例がいくつか見つかりません。

def isGrafting(a):
  for i in xrange(1, int(ceil(log10(a))) + 2):
    if a == floor((sqrt(a) * 10**(i-1)) % 10**int(ceil(log10(a)))):
      return 1

a = 0
while(1):
  if (isGrafting(a)):
    print "%d %.15f" % (a, sqrt(a))
  a += 1

このコードには、少なくとも 1 つの既知の移植番号がありません。9999999998 => 99999.99998999999999949999999994999999999374999999912... を掛けた後、余分な精度が落ちるよう10**5です。

>>> a = 9999999998
>>> sqrt(a)
99999.99999
>>> a == floor((sqrt(a) * 10**(5)) % 10**int(ceil(log10(a))))
False
>>> floor((sqrt(a) * 10**(5)) % 10**int(ceil(log10(a))))
9999999999.0
>>> print "%.15f" % sqrt(a)
99999.999989999996615
>>> print "%.15f" % (sqrt(a) * 10**5)
9999999999.000000000000000

そこで、短い C++ プログラムを作成して、CPU が浮動小数点数を切り捨てたのか、何らかの方法で python を切り捨てたのかを確認しました。

#include <cstdio>
#include <cmath>
#include <stdint.h>

int main()
{
  uint64_t a = 9999999998;
  printf("%ld %.15f %.15f %.15f %.15f\n", a, sqrt((double)a), sqrt((double)a)*1e4, sqrt((double)a)*1e5, sqrt((double)a)*1e6);
  a = 999999999998;
  printf("%ld %.15f %.15f %.15f %.15f\n", a, sqrt((double)a), sqrt((double)a)*1e5, sqrt((double)a)*1e6, sqrt((double)a)*1e7);
  a = 99999999999998;
  printf("%ld %.15f %.15f %.15f %.15f\n", a, sqrt((double)a), sqrt((double)a)*1e6, sqrt((double)a)*1e7, sqrt((double)a)*1e8);
  return 0;
}

どの出力:

9999999998 99999.999989999996615 999999999.899999976158142 9999999999.000000000000000 99999999990.000000000000000
999999999998 999999.999998999992386 99999999999.899993896484375 999999999999.000000000000000 9999999999990.000000000000000
99999999999998 9999999.999999899417162 9999999999999.900390625000000 99999999999999.000000000000000 999999999999990.000000000000000

そのため、浮動小数点精度の制限に対して激しく実行しているように見え、CPU は残りの差が浮動小数点エラーであると考えているため、残りのビットを切り落としています。Python でこれを回避する方法はありますか? または、C に移行して GMP などを使用する必要がありますか?

4

5 に答える 5

57

標準ライブラリでは、decimalモジュールが探しているものである可能性があります。また、mpmathが非常に役立つことがわかりました。ドキュメントには多くの優れた例も含まれています (残念ながら、私のオフィスのコンピューターにはインストールされていませんmpmath。それ以外の場合は、いくつかの例を確認して投稿します)。

decimalただし、モジュールに関する注意点が 1 つあります。モジュールには、単純な数学演算用の組み込み関数がいくつか含まれています (例: sqrt)。ただし、これらの関数の結果は、または他のモジュールの対応する関数とmathより高い精度で常に一致するとは限りません (より正確な場合もあります)。例えば、

from decimal import *
import math

getcontext().prec = 30
num = Decimal(1) / Decimal(7)

print("   math.sqrt: {0}".format(Decimal(math.sqrt(num))))
print("decimal.sqrt: {0}".format(num.sqrt()))

Python 3.2.3 では、これは最初の 2 行を出力します。

   math.sqrt: 0.37796447300922719758631274089566431939601898193359375
decimal.sqrt: 0.377964473009227227214516536234
actual value: 0.3779644730092272272145165362341800608157513118689214

述べたように、これは期待どおりではなく、精度が高いほど結果が一致しないことがわかります。この例では、実際の値decimalとより厳密に一致するため、モジュールの精度がより高いことに注意してください。

于 2012-07-17T13:05:25.493 に答える
11

この特定の問題についてdecimalは、10 進数をタプルとして格納するため、 が最適な方法です。

>>> a = decimal.Decimal(9999999998)
>>> a.as_tuple()
DecimalTuple(sign=0, digits=(9, 9, 9, 9, 9, 9, 9, 9, 9, 8), exponent=0)

10 進表記で最も自然に表現されるプロパティを探しているので、2 進表記を使用するのは少しばかげています。リンク先のウィキペディアのページには、「グラフト数字」が始まる前に表示される「グラフトしない数字」の数が示されていなかったため、次のように指定できます。

>>> def isGrafting(dec, max_offset=5):
...     dec_digits = dec.as_tuple().digits
...     sqrt_digits = dec.sqrt().as_tuple().digits
...     windows = [sqrt_digits[o:o + len(dec_digits)] for o in range(max_offset)]
...     return dec_digits in windows
... 
>>> isGrafting(decimal.Decimal(9999999998))
True
>>> isGrafting(decimal.Decimal(77))
True

2 進数表現と 10 進数表現の間の変換Decimal.sqrt()のため、 の結果は の結果よりも、少なくともこれについてはより正確になる可能性が高いと思います。math.sqrt()たとえば、次のことを考慮してください。

>>> num = decimal.Decimal(1) / decimal.Decimal(7)
>>> decimal.Decimal(math.sqrt(num) ** 2) * 7
Decimal('0.9999999999999997501998194593')
>>> decimal.Decimal(num.sqrt() ** 2) * 7
Decimal('1.000000000000000000000000000')
于 2012-07-17T13:24:35.977 に答える
9

浮動小数点の代わりにDecimalを試すことができます。

于 2012-07-17T12:57:30.350 に答える
6

Python には組み込みの任意精度浮動小数点数はありませんが、GMP を使用するサードパーティの Python パッケージがあります: gmpyおよびPyGMP

于 2012-07-17T13:02:49.457 に答える
2

を使用decimalします (より明確な例を次に示します):

>>> 2.3-2.2
0.09999999999999964
>>> from decimal import Decimal
>>> Decimal('2.3')-Decimal('2.2')
Decimal('0.1')
>>> float(Decimal('2.3')-Decimal('2.2'))
0.1
>>> 
于 2018-09-18T01:29:41.540 に答える