具有mpz/mpfr值的numpy数组

Question

我想要一个带有mpz/mpfr值的numpy数组。因为我的代码：

import numpy as np
import gmpy2
A=np.ones((5,5));
print A/gmpy2.mpfr(1);

生成：

RuntimeWarning: invalid value encountered in divide
  print A/gmpy2.mpfr(1);
[[1.0 1.0 1.0 1.0 1.0]
 [1.0 1.0 1.0 1.0 1.0]
 [1.0 1.0 1.0 1.0 1.0]
 [1.0 1.0 1.0 1.0 1.0]
 [1.0 1.0 1.0 1.0 1.0]]

据我所知，将gmpy mpfr转换为numpy float64是不可能的。那么，我如何才能得到一个带有mpfr值的numpy数组呢？

谢谢。

Answer 1

您需要使用dtype=object创建数组，然后可以在数组中使用任何python类型。我没有安装gmpy2，但下面的示例应该说明它是如何工作的：

In [3]: a = np.ones((5, 5), dtype=object)

In [5]: import fractions

In [6]: a *= fractions.Fraction(3, 4)

In [7]: a
Out[7]: 
array([[3/4, 3/4, 3/4, 3/4, 3/4],
       [3/4, 3/4, 3/4, 3/4, 3/4],
       [3/4, 3/4, 3/4, 3/4, 3/4],
       [3/4, 3/4, 3/4, 3/4, 3/4],
       [3/4, 3/4, 3/4, 3/4, 3/4]], dtype=object)

拥有dtype=object的numpy数组可能会有一点误导，因为让标准dtype操作变得超级快的强大numpy机制现在由默认对象的python操作符负责，这意味着速度将不再存在：

In [12]: b = np.ones((5, 5)) * 0.75

In [13]: %timeit np.sum(a)
1000 loops, best of 3: 1.25 ms per loop

In [14]: %timeit np.sum(b)
10000 loops, best of 3: 23.9 us per loop

Answer 2

我相信这是两个库中的一个错误。我也相信它是被修复的。

输入：

import sys
import numpy as np
import gmpy2

print(sys.version)
print(np.__version__)
print(gmpy2.version)

A=np.ones((5,5));
print(A/gmpy2.mpfr(1))

输出：

3.4.2 (v3.4.2:ab2c023a9432, Oct  6 2014, 22:15:05) [MSC v.1600 32 bit (Intel)]
1.9.1
2.0.5
[[mpfr('1.0') mpfr('1.0') mpfr('1.0') mpfr('1.0') mpfr('1.0')]
 [mpfr('1.0') mpfr('1.0') mpfr('1.0') mpfr('1.0') mpfr('1.0')]
 [mpfr('1.0') mpfr('1.0') mpfr('1.0') mpfr('1.0') mpfr('1.0')]
 [mpfr('1.0') mpfr('1.0') mpfr('1.0') mpfr('1.0') mpfr('1.0')]
 [mpfr('1.0') mpfr('1.0') mpfr('1.0') mpfr('1.0') mpfr('1.0')]]

当遇到未知类型时，Numpy没有正确地说明该怎么做，或者gmpy2没有指定如何除以某物(__rdiv__)。

没有必要指定ndarray的dtype，除非您打算覆盖它的元素。像乘法这样的操作会产生一个新的ndarray，Numpy会找出要使用的dtype。