我有一个形状为M*3的大型数组A,它的每一行的元素都是从0到N-1的唯一非负整数.实际上,在我的有限元分析中,每一行对应一个三角形.

例如,M=4,N=5,矩阵A如下所示

array([[0, 1, 2],
       [0, 2, 3],
       [1, 2, 4],
       [3, 2, 4]])

现在,我需要构造另一个大小为M*N的数组B,以便

B[m,n] = 1 if n is in A[m], or else 0 

上述示例性A的对应B将是

1 1 1 0 0
1 0 1 1 0
0 1 1 0 1
0 0 1 1 1

基于循环的代码将是

B = np.zeros((M,N))
for m in range(M):
  for n in B[m]:
    B[m,n]=1

但是由于我有很大的M和N(每个都有10^6的刻度),我如何使用好的块索引技术来加速这个过程?此外,我觉得也需要稀疏矩阵技术,因为1字节的M*N数据大约是10**12,即1000G.

总的来说,我觉得使用NumPy的矢量化技术,如索引和广播,看起来更像是一种临时的、容易出错的活动,依赖于相当多的街头智慧(如果你愿意,也可以称之为艺术).有没有什么编程语言可以系统地将基于循环的代码转换为高性能的矢量化版本?

推荐答案

您可以根据您的数据直接创建稀疏CSR矩阵

正如您在问题中已经提到的,由uint8值组成的密集矩阵需要1 TB.通过使用稀疏矩阵,可以将其减少到大约.19 MB,如下面的示例所示.

Creating Inputs with relevant size

这应该包括在问题中,因为它给出了矩阵的稀疏性的提示.

from scipy import sparse
import numpy as np

M=int(1e6)
N=int(1e6)

A=np.random.randint(low=0,high=N,size=(M,3))

Creating a sparse csr-matrix

看一看这scipy-doc篇文章,或者想要了解一下一般性的概述,wiki篇文章也可能会很有用.

#array of ones, the same size of non-zero values (3 MB if uint8)
data   =np.ones(A.size,dtype=np.uint8)

#you already have the indices, they are expected as an 1D-array (12 MB)
indices=A.reshape(-1)

#every A.shape[1], a new row beginns (4 MB)
indptr =np.arange(0,A.size+1,A.shape[1])

B=sparse.csr_matrix((data,indices,indptr),shape=(M, N))

Python相关问答推荐

过滤绕轴旋转的螺旋桨

Python -Polars库中的滚动索引?

返回nxon矩阵的diag元素,而不使用for循环

如何在BeautifulSoup中链接Find()方法并处理无?

通过优化空间在Python中的饼图中添加标签

在Python中处理大量CSV文件中的数据

难以在Manim中正确定位对象

Pandas - groupby字符串字段并按时间范围 Select

如何获取TFIDF Transformer中的值?

如何在Python数据框架中加速序列的符号化

如何使Matplotlib标题以图形为中心,而图例框则以图形为中心

如何指定列数据类型

在Python 3中,如何让客户端打开一个套接字到服务器,发送一行JSON编码的数据,读回一行JSON编码的数据,然后继续?

如何在Pyplot表中舍入值

将标签移动到matplotlib饼图中楔形块的开始处

Python—在嵌套列表中添加相同索引的元素,然后计算平均值

提取最内层嵌套链接

每次查询的流通股数量

Pandas数据框上的滚动平均值,其中平均值的中心基于另一数据框的时间

Pandas:计数器的滚动和,复位