当使用无效的numpy数组时返回西皮.利纳格.布拉斯.sgemm公司

Question

我在计算A•A^T：

# These are my dummy values for testing
A = np.ones((150000,265),dtype=np.float32, order='F')
A_T = np.ones((265, 150000),dtype=np.float32, order='F')

out = scipy.linalg.blas.sgemm(alpha=1.0, a=A, b=A_T)

两分钟后：

^{pr2}$

注意到0了吗？我迷路了。。。我尝试过使用64位浮点并得到相同的输出。 从35468开始，数组为零。在

In [39]: out[0,35468]
Out[39]: 0.0

In [9]: scipy.__version__
Out[9]: '0.12.1'

更新/编辑：

我相当肯定美国运输部正在调用*gemm方法本身。在

In [1]: A = np.ones((150000,265), dtype=np.float32, order='F')

In [2]: A_T = np.ones((265, 150000),dtype=np.float32, order='F')

In [3]: out = A.dot(A_T)

In [4]: out.shape
Out[4]: (150000, 150000)

In [5]: out
Out[5]: 
array([[ 265.,  265.,  265., ...,  265.,  265.,  265.],
   [ 265.,  265.,  265., ...,  265.,  265.,  265.],
   [ 265.,  265.,  265., ...,  265.,  265.,  265.],
   ..., 
   [   0.,    0.,    0., ...,    0.,    0.,    0.],
   [   0.,    0.,    0., ...,    0.,    0.,    0.],
   [   0.,    0.,    0., ...,    0.,    0.,    0.]], dtype=float32)

Answer 1

在我的机器上，您的示例以Memory Error死亡，与{}相同。在

在减小矩阵尺寸时，它产生一个正确的结果。在

要进一步调试它，请尝试直接使用Fortran BLAS调用，而不使用python。如果通过，请考虑filing a bug in the scipy tracker。在

另外，从scipy0.13以后，您可以使用syrk。在