归一化向量在 Numpy 中产生 nan答案

【问题标题】：Normalizing vector produces nan in Numpy归一化向量在 Numpy 中产生 nan
【发布时间】：2015-07-29 21:53:56
【问题描述】：

我从 scipy/numpy 得到一些奇怪的行为，我怀疑这是一个错误，但有人可能更了解？我有一对长数组，为了调试目的，我将它们分成长度为 2-4 的帧。我想标准化每对帧并取点积。执行此操作的代码（带有一些调试输出）是：

   tf = numpy.copy(t_frame) / norm(t_frame)
   pf = numpy.copy(p_frame) / norm(p_frame)
   print "OPF:"
   print p_frame
   print "PF: "
   print pf
   print "TF norm is: " + str(norm(tf))
   print "PF norm is: " + str(norm(pf))
   print numpy.dot(tf, pf)
   return numpy.dot(tf, pf)

这符合我一段时间的预期（特别是 tf 和 pf 的标准为 1），但随后我开始看到这样的行：

OPF：

[-91 -119 -137 -132]

PF：

[楠楠楠楠]

什么？这可以在新的 Python 窗口中正常化：

>>> p = [ -91, -119, -137, -132] 
>>> p / norm(p)
array([-0.37580532, -0.49143773, -0.56577285, -0.54512421])

为了它的价值，我尝试了numpy.linalg.norm、scipy.linalg.norm，并定义了一个函数来返回点积的平方根。

有什么想法吗？

更新：感谢您的建议！我尝试将 dtype 切换为 float128，但遗憾的是得到了类似的行为。我实际上倾向于认为这是 Python 中的一个错误，而不是此时的 numpy：

如果这是一个直截了当的溢出问题，似乎我会得到与给定列表一致的问题。但是，如果我在新的 python 会话中执行它，规范计算得很好。
我尝试自己滚动：
```
def norm(v):
   return (  sum(numpy.array(v)*numpy.array(v)))**(0.5)
```
这仅使用 numpy 来表示数组。我仍然遇到同样的问题，但 later 在数据集中（并且没有运行时警告）。它正在执行大约 37000 次这样的计算。
我实际上是在两个帧上计算标准，一个 t_frame 和一个 p_frame。当且仅当另一个计算发生时，一个的计算才会阻塞。

综上所述，我认为 Python (2.7.9) 内部某处存在一些奇怪的缓冲区溢出？？？我最终也需要这些计算快速；所以我正在考虑切换到 Cython 进行计算。

更新 2：我真的尝试过自己动手：

def norm(v):
  sum = float(0)
  for i in range(len(v)):
    sum += v[i]**2
  return sum**(0.5)

问题就消失了。所以我猜它是 numpy 中的一个错误（Gentoo Linux 上的 1.9.0）。

【问题讨论】：

标签： python numpy scipy

【解决方案1】：

看起来这是 numpy 中的一个错误。如果数组的数据类型是np.int16，我可以重现问题：

In [1]: np.__version__
Out[1]: '1.9.2'

In [2]: x = np.array([ -91, -119, -137, -132], dtype=np.int16)

In [3]: x
Out[3]: array([ -91, -119, -137, -132], dtype=int16)

In [4]: np.linalg.norm(x)
/Users/warren/anaconda/lib/python2.7/site-packages/numpy/linalg/linalg.py:2061: RuntimeWarning: invalid value encountered in sqrt
  return sqrt(sqnorm)
Out[4]: nan

这个问题也出现在开发版numpy的master分支上。我在这里创建了一个问题：https://github.com/numpy/numpy/issues/6128

如果p_frame 实际上是一个 16 位整数数组，一个简单的解决方法是：

x = np.asarray(p_frame, dtype=np.float64)
pf = x / norm(x)

【讨论】：

【解决方案2】：

按照 Warren 的链接之一，我收到以下警告：

In [1016]: np.linalg.norm(100000*np.ones(2).astype('int16'))
/usr/local/lib/python2.7/site-packages/numpy/linalg/linalg.py:2051: RuntimeWarning: invalid value encountered in sqrt
  return sqrt(add.reduce((x.conj() * x).real, axis=None))

对于这个x2，内部表达式是否定的——小dtype中溢出的结果。

In [1040]: x2=100000*np.ones(2).astype('int16')
In [1041]: np.add.reduce((x2.conj()*x2).real,axis=None)
Out[1041]: -1474836480

与x1类似：

In [1042]: x1
Out[1042]: array([ -9100, -11900, -13700, -13200], dtype=int16)
In [1043]: np.add.reduce((x1.conj()*x1).real,axis=None)
Out[1043]: -66128

如果“点”的总和对于dtype 来说太大，它可能是负数，当通过sqrt 时产生nan。

（我在linux下使用的是1.8.2和1.9.0）。

【讨论】：