为什么np.dot不精确？（n维数组）

Question

假设我们采用np.dot两个'float32'2D数组：

res = np.dot(a, b)   # see CASE 1
print(list(res[0]))  # list shows more digits

[-0.90448684, -1.1708503, 0.907136, 3.5594249, 1.1374011, -1.3826287]

数字。除了它们可以更改：

---

案例1：切片a

np.random.seed(1)
a = np.random.randn(9, 6).astype('float32')
b = np.random.randn(6, 6).astype('float32')

for i in range(1, len(a)):
    print(list(np.dot(a[:i], b)[0])) # full shape: (i, 6)

[-0.9044868,  -1.1708502, 0.90713596, 3.5594249, 1.1374012, -1.3826287]
[-0.90448684, -1.1708503, 0.9071359,  3.5594249, 1.1374011, -1.3826288]
[-0.90448684, -1.1708503, 0.9071359,  3.5594249, 1.1374011, -1.3826288]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]

即使打印的切片从完全相同的数字相乘得出的结果也有所不同。

---

案例2：展平a，取的一维版本b，然后切片a：

[-0.9044868,  -1.1708502, 0.90713596, 3.5594249, 1.1374012, -1.3826287]
[-0.90448684, -1.1708503, 0.9071359,  3.5594249, 1.1374011, -1.3826288]
[-0.90448684, -1.1708503, 0.9071359,  3.5594249, 1.1374011, -1.3826288]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]
[-0.90448684, -1.1708503, 0.907136,   3.5594249, 1.1374011, -1.3826287]

np.random.seed(1)
a = np.random.randn(9, 6).astype('float32')
b = np.random.randn(1, 6).astype('float32')

for i in range(1, len(a)):
    a_flat = np.expand_dims(a[:i].flatten(), -1) # keep 2D
    print(list(np.dot(a_flat, b)[0])) # full shape: (i*6, 6)

---

案例3：更强的控制力；将所有不涉及的整数设置为零：添加a[1:] = 0到CASE 1代码中。结果：差异仍然存在。

---

案例4：检查索引以外的内容[0]; 与for一样[0]，结果从其创建点开始就稳定了固定数量的数组扩展。输出量

[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]
[-0.3393164, 0.9528787, 1.3627989, 1.5124314, 0.46389243, 1.437775]

---

因此，对于2D * 2D情况，结果有所不同-但对于1D * 1D却是一致的。从我的一些读物来看，这似乎源于使用简单加法的1D-1D，而2D-2D使用了“精确”的性能提升加法，但精度可能较低（例如，成对加法则相反）。但是，我无法理解为什么一旦情况1 a越过设定的“阈值”，差异就会消失。大a和b，后来这个门槛似乎撒谎，但它始终存在。

所有人都说：为什么np.dotND-ND阵列不精确（且不一致）？相关的Git

---

附加信息：

• 环境：Win-10操作系统，Python 3.7.4，Spyder 3.3.6 IDE，Anaconda 3.0 2019/10
• CPU：i7-7700HQ为2.8 GHz
• Numpy v1.16.5

罪魁祸首可能是库：NumPy的MKL -也BLASS库; 感谢Bi Rico的注意

---

压力测试代码：如前所述，与较大的阵列相比，频率差异会加剧；如果上面无法再现，则下面应该是可再现的（如果不能再现，请尝试更大的暗色）。我的输出

np.random.seed(1)
a = np.random.randn(9, 6).astype('float32')
b = np.random.randn(6, 6).astype('float32')

for j in range(len(a) - 2):
    for i in range(1, len(a)):
        res = np.dot(a[:i], b)
        try:    print(list(res[j]))
        except: pass
    print()

---

问题严重性：显示的差异“很小”，但在神经网络上运行时不再如此，数十亿个数字在几秒钟内成倍增加，而在整个运行过程中则为数万亿个。每个线程报告的模型精度相差百分之十。

下面是将模型馈入模型后得到的数组的gif a[0]，w / len(a)==1vs len(a)==32：

---

根据Paul的测试，并感谢其他平台的结果：

案例1转载（部分）：

• Google Colab VM-英特尔至强2.3 G-Hz-Jupyter-Python 3.6.8
• Win-10 Pro Docker桌面-英特尔i7-8700K-jupyter / scipy-notebook-Python 3.7.3
• Ubuntu 18.04.2 LTS + Docker-AMD FX-8150-jupyter / scipy-notebook-Python 3.7.3

注意：这些产生的错误比上面显示的要低得多；第一行中的两个条目与其他行中的相应条目的最低有效位相差1。

情况1未转载：

• Ubuntu 18.04.3 LTS-Intel i7-8700K-IPython 5.5.0-Python 2.7.15+和3.6.8（2个测试）
• Ubuntu 18.04.3 LTS-英特尔i5-3320M-IPython 5.5.0-Python 2.7.15+
• Ubuntu 18.04.2 LTS-AMD FX-8150-IPython 5.5.0-Python 2.7.15rc1

注意事项：

• 该链接 Colab笔记本和jupyter环境中表现出远较小的差异（只为前两行）比我的系统上观察到。同样，案例2从未（至今）显示出不精确性。
• 在这个非常有限的示例中，当前（Dockerized）的Jupyter环境比IPython环境更容易受到攻击。
• np.show_config()发布时间太长，但总而言之：IPython环境是基于BLAS / LAPACK的；Colab基于OpenBLAS。在IPython Linux环境中，BLAS库是系统安装的-在Jupyter和Colab中，它们来自/ opt / conda / lib

---

更新：可接受的答案是准确的，但范围广泛且不完整。对于任何可以在代码级别解释行为的人来说，这个问题仍然悬而未决-即，所使用的精确算法np.dot，以及如何解释上述结果中观察到的“一致不一致”（另请参见注释）。这是我无法理解的一些直接实现：sdot.c - arraytypes.c.src

Answer 1

这看起来是不可避免的数值不精确性。如所解释这里，NumPy的使用对于矩阵乘法高度优化，仔细调谐BLAS方法。这意味着可能要乘以两个矩阵的运算顺序（总和与乘积）会随着矩阵大小的变化而变化。

为了更加清晰，我们知道，从数学上讲，可以将所得矩阵的每个元素计算为两个向量（数字的等长序列）的点积。但这不是 NumPy计算所得矩阵元素的方式。实际上，有更有效但复杂的算法，例如Strassen算法，无需直接计算行列点积即可获得相同的结果。

使用此类算法时，即使将所得矩阵C = AB的元素C ij在数学上定义为A 的第i行与B 的第j列的点积，如果您将矩阵A2与相同的第i行作为阿与基质B2具有相同的第j列乙，元件C2 IJ将实际计算以下操作中的一个不同的序列（即依赖于整个A2和B2 矩阵），可能会导致不同的数值误差。

这就是为什么即使在数学上C ij = C2 ij（如在案例1中一样），计算中算法遵循的不同操作顺序（由于矩阵大小的变化）也会导致不同的数值误差。数值误差还解释了根据环境以及在某些情况下对于某些环境可能不存在数值误差这一事实而略有不同的结果。