Numpy高效大矩阵乘法

Question

要在磁盘上存储大矩阵,我使用numpy.memmap.

以下是测试大矩阵乘法的示例代码:

import numpy as np
import time

rows= 10000 # it can be large for example 1kk
cols= 1000

#create some data in memory
data = np.arange(rows*cols, dtype='float32') 
data.resize((rows,cols))

#create file on disk
fp0 = np.memmap('C:/data_0', dtype='float32', mode='w+', shape=(rows,cols))
fp1 = np.memmap('C:/data_1', dtype='float32', mode='w+', shape=(rows,cols))

fp0[:]=data[:]
fp1[:]=data[:]

#matrix transpose test
tr = np.memmap('C:/data_tr', dtype='float32', mode='w+', shape=(cols,rows))
tr= np.transpose(fp1)  #memory consumption?
print fp1.shape
print tr.shape

res = np.memmap('C:/data_res', dtype='float32', mode='w+', shape=(rows,rows))
t0 = time.time()
# redifinition ? res= np.dot(fp0,tr) #takes 342 seconds on my machine, if I multiplicate matrices in RAM it takes 345 seconds (I thinks it's a strange result)
res[:]= np.dot(fp0,tr) # assignment ?
print res.shape
print (time.time() - t0)

所以我的问题是:

1. 如何限制使用此过程的应用程序的内存消耗,例如100Mb(或1Gb或其他)的某些值.还有我不明白如何估计程序的内存消耗(我认为内存只在"数据时分配" "创建变量,但是当我们使用memmap文件时会使用多少内存？)
2. 也许存在一些存储在磁盘上的大矩阵乘法的最优解？例如,数据可能没有最佳地存储在磁盘上或者从磁盘中获取,没有正确处理,并且点产品只使用一个核心.也许我应该使用像PyTables这样的东西？

我也对使用有限内存使用的线性方程组(SVD和其他)的算法感兴趣.也许这种算法称为核外或迭代,我认为有一些类比,如硬盘< - > ram,gpu ram < - > cpu ram,cpu ram < - > cpu cache.

同时在这里我找到了关于PyTables矩阵乘法一些信息.

我也在R中发现了这个,但我需要它用于Python或Matlab.

Answer 1

Dask.array使用阻塞算法和任务调度为大型磁盘阵列提供了一个numpy接口.它可以很容易地进行核外矩阵乘法和其他简单的numpy操作.

阻塞的线性代数更难,你可能想看看关于这个主题的一些学术着作.Dask确实支持高和瘦的矩阵上的QR和SVD分解.

无论是大型数组,你真的想要阻塞算法,而不是天真的遍历,这会以不愉快的方式击中磁盘.

Answer 2

考虑使用NumExpr进行处理：https://github.com/pydata/numexpr

\n\n

...在内部，NumExpr 采用自己的矢量化虚拟机，该虚拟机是围绕分块读取策略设计的，以便有效地操作内存中最佳大小的数据块。如果调整得当，它可以轻松击败 na\xc3\xafve NumPy 操作。

\n\n

NumExpr 可能涵盖了问题细分中的#2。如果您通过使用可流式二进制格式解决#1，则可以在加载数据文件 \xe2\x80\x93 时采用分块读取方法，如下所示：

\n\n

    with open(\'path/to/your-data.bin\', \'rb\') as binary:\n        while True:\n            chunk = binary.read(4096) # or what have you\n            if not chunk:\n                break\n

\n\n

如果这对您来说太低级，我建议您查看HDF5库和格式：http://www.h5py.org \xe2\x80\x93 it\xe2\x80\x99s 二进制文件的最佳解决方案据我所知，基于 NumPy 的结构的序列化。该h5py模块支持压缩、分块读取、dtypes、元数据\xe2\x80\xa6 等等。

\n\n

祝你好运！

\n