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

此代码至少缺少一个已知的嫁接编号。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不过，关于该模块的一个警告。该模块包含几个用于简单数学运算的内置函数（例如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 中，这会输出前两行

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

如前所述，这并不是您所期望的，您可以看到精度越高，结果匹配越少。请注意，该decimal模块在此示例中确实具有更高的准确性，因为它与实际值更接近。

于 2012-07-17T13:05:25.493 回答
11

对于这个特殊的问题，这decimal是一个很好的方法，因为它将十进制数字存储为元组！

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

由于您正在寻找最自然地以十进制表示的属性，因此使用二进制表示有点愚蠢。您链接到的维基百科页面没有说明在“嫁接数字”开始之前可能出现多少个“非嫁接数字”，因此您可以指定：

>>> 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

我认为由于二进制表示和十进制表示之间的转换，结果很有可能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 的第 3 方 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 回答