如何在python中有效地计算（稀疏）位矩阵的矩阵乘积

Posted 2023-03-12

【中文标题】如何在python中有效地计算（稀疏）位矩阵的矩阵乘积

【英文标题】：How to calculate matrix product on (sparse) bit matrix efficiently in python

【发布时间】：2016-10-01 03:52:50

【问题描述】：

位数组之间的普通矩阵乘积：

           1, 0, 0
Matrix A = 1, 1, 1
           0, 0, 1

                 1, 1, 0
Transpose of A = 0, 1, 0
                 0, 1, 1

C = Matrix A times (Transpose of Matrix A)

    1, 1, 0
C = 1, 3, 1
    0, 1, 1

A 是 1 和 0 的位数组。实际矩阵 A 很大，大约有 0.25% 的 1 项和 99.75% 的零项。

C 是一个整数数组。

如何在不占用大量内存的情况下快速计算？

目前，我正在使用 scipy 的稀疏矩阵乘法例程来处理 python 中浮点 1.0 和 0.0 的压缩稀疏行矩阵。我也在尝试直接调用 mkl 库中的 c 函数来减少内存使用。

【问题讨论】：

多大？ 1000 x 1000 矩阵需要多长时间？

计算一个 750,000 x b 矩阵与其转置的乘积。 b 至少为 10,000。 b 越大越好。但是每个 cpu 核心的内存是有限的。

对不起，我的意思是一个形状为 b x 750,000 的矩阵 A。所以 A 的 A 次转置应该是 b x b 的形状。

@rxu 给出 750000 和 0.25%，你的 C 有 99% 的非零元素。您可能希望以密集格式存储它。

这取决于 b 行是否倾向于在同一列有 1。

【参考方案1】：

MKL 等现有性能库始终使用 float/double 作为数据类型。与将 A 转换为浮点 CSR 然后调用 .dot() 或一些 MKL 例程相比，您可能会发现编写自己的 bit-mat-mul 代码更快。你甚至不需要乘法运算。它只是计算位数。

编辑

在了解您的问题背景后，我会建议以下程序。

将数组 A 转换为 CSR 格式并仅存储 col 索引和行 ptr； 对于 A 的每一行 i 和 j 行，计算公共 col 索引的数量并将结果存储在 C(i, j) 中，其中 i

稠密矩阵 C 就是你想要的。

给定你的 A 的大小 (b x 750,000) 和密度 (0.25%)，C 的密度是 99.1%；每个 col 索引的平均长度为 1875。

所以你的问题变成了计算 2 个 1875-D 向量的公共元素的数量 b*(b+1)/2 次。

for 循环的速度似乎是唯一剩下的问题。

【讨论】：

@rxu 之前的问题暗示数据可能存储为gmpy2.mpz 类型。如果正确，gmpy2.hamdist(row & col, 0) 应该很快。

让我试试。我担心开销 python 循环和从 python 调用 c++ 会使事情变慢......实际上我写了一个版本来使用 python 的集合和交集一次。对于稀疏矩阵，那个比 scipy 的点要慢得多。即使构建稀疏矩阵需要时间。组合时间仍然比设置版本快得多。在程序的前一部分中，我使用 mpz 和 python 循环非常快速地计算了logical_and。效果很好，而 scipy 做不到那么快，而且会占用大量内存。

在密集位数组的程序工作的第一部分。现在我必须在程序的第二部分使用稀疏数组。

set 版本也占用大量内存。这使得 b 很小。

另一件事是我必须在计算产品之前将该 b x 750,000 矩阵压缩到内存中。在我的例子中，稀疏矩阵在内存效率上大大优于密集矩阵。