如何避免在scipy稀疏矩阵中插入不必要的零

Question

我在玩弄一个由线性系统建模的问题，它可以写成一个方块三对角矩阵。这些块的尺寸为b=4n+8，整个矩阵的尺寸为Nb；N可以任意大（当然是合理的），而N保持相当小（通常小于10）

块本身是稀疏的，第一条对角线仅为单位矩阵，第二条对角线每个块仅具有n+1个非零列（因此3n+7列零）。这些列是连续的，要么是零，要么是非零，要么是相反

在内存中构建所有这些块会产生一个3N-2 x b x b数组，该数组可以使用scipy.sparse.bsr_矩阵转换为稀疏矩阵，然后转换为CSR格式并修剪多余的零。它工作得很好，但我宁愿跳过这个临时的大型稀疏数组（对于N=1e4，N=5，每个相关条目都是5.6个零！）

• 我看了一下scipy.sparse.dok_矩阵，建议用于切片和增量构建。创建我的条目适合于一个整洁的循环，但是这个过程比使用bsr_矩阵和我不必要的密集数组要长约10倍，这将对未来的用例有害
• bsr_矩阵似乎不能直接用于scipy稀疏矩阵作为输入
• 使用bsr_矩阵而不使用包括对角块，然后添加稀疏眼大大减少了零的数量（在我的测试配置中每个相关条目3.5个）并且与原始解决方案相比，过程加快了三分之一。得分

关于我可以做些什么来进一步减少这个矩阵的原始印记，有什么线索吗？显而易见的目标是给我更多的自由选择N

编辑

通过分别构造三个块对角线，我设法将事情做了一些改进。通过这样做，我的第二条对角线需要更少的填充（n+3而不是3n+7；每个相关条目1.3个零），将我的原始块划分为两个垂直块（一个充满零）和我一次只需要一条对角线，在此基础上将内存成本减半。主对角线仍然是用眼法构造的。锦上添花：与我的第三个要点相比，速度提高了25%，这可能是因为分离两条第二对角线节省了使用bsr_矩阵之前所需的一些阵列重塑操作。与原始方法相比，对于我的（N，N）=（1e4，5）测试用例，在修剪之前比较矩阵时节省了约2000万个零。每个128位，已经是不错的增益了

我现在能想象到的唯一可能的改进是分别构建这些对角线，没有任何填充，然后插入零列（可能通过具有标识块矩阵的乘积），最后将所有内容添加到一起。 我还读了一些关于使用dict更新空dok_矩阵的内容，但在我的例子中，我认为我需要扩展索引列表，使用它们的笛卡尔积来构造键，我的块的每个元素都需要是一个单独的值，因为显然不能使用切片作为字典键

Answer 1

我最终实现了我在最后一段中提出的解决方案。 对于每个第二对角线，我构造一个没有任何填充的块稀疏矩阵，然后通过一个块为单位的块矩阵的右手边积将其转换为一个适当形状的矩阵。我需要在这里存储零以使用bsr_矩阵（我首先尝试了scipy.sparse.block_diag方法，但速度非常慢），但与我的半填充解决方案相比，它们的数量更少：（4n+7）（n+1）vs（4n+8）（n+3）；它们可以用8位而不是128位来表示。执行时间增加了约40%，但我可以接受（与第一个解决方案相比，它仍然减少了20%）

我可能遗漏了一些东西，但现在我对这个解决方案非常满意

编辑

在影响乘积之前修剪RHS矩阵的零点时，与之前最有效的解决方案相比，执行时间减少30%；结局好，一切都好