如何重复使用 block_diag答案

【问题标题】：How to use block_diag repeatedly如何重复使用 block_diag
【发布时间】：2016-03-11 08:53:14
【问题描述】：

我有一个相当简单的问题，但仍然无法解决。我想要一个块对角线 n^2*n^2 矩阵。这些块是稀疏的 n*n 矩阵，只有对角线，首先是对角线，然后是对角线。对于n=4 的简单情况，这很容易做到

datanew = ones((5,n1))
datanew[2] = -2*datanew[2]
diagsn = [-4,-1,0,1,4]
DD2 = sparse.spdiags(datanew,diagsn,n,n)
new = sparse.block_diag([DD2,DD2,DD2,DD2])

由于这只对小 n 有用，有没有更好的方法来使用 block_diag？考虑 n -> 1000

【问题讨论】：

我稍微清理了你的显示器。您需要更具体地说明为什么这不适用于较大的n。（n1 是什么？）block_diag 和 bmat 是开放的 Python 代码，因此您可以研究它们，如果需要简化操作以适应您的情况。
bmat 最终将您的 DD2 转换为 coo 格式，将他们的 data、row、col 加入到 3 个大数组中，并从那个。
重点是我不想写block_diag(DD2,....... (DD2,:) 直到 n) ?
一种粗略的方法是将DD2 附加到列表中一千次，然后将该列表传递给block_diag。
这可以在 block_diag 中使用吗？我不确定它是否会接受列表

标签： python-2.7 scipy sparse-matrix

【解决方案1】：

构造一长串DD2 矩阵的简单方法是使用列表推导：

In [128]: sparse.block_diag([DD2 for _ in range(20)]).A
Out[128]: 
array([[-2,  1,  0, ...,  0,  0,  0],
       [ 1, -2,  1, ...,  0,  0,  0],
       [ 0,  1, -2, ...,  0,  0,  0],
       ..., 
       [ 0,  0,  0, ..., -2,  1,  0],
       [ 0,  0,  0, ...,  1, -2,  1],
       [ 0,  0,  0, ...,  0,  1, -2]])

In [129]: _.shape
Out[129]: (80, 80)

至少在我的版本中，block_diag 想要一个数组列表，而不是 *args：

In [133]: sparse.block_diag(DD2,DD2,DD2,DD2)
...
TypeError: block_diag() takes at most 3 arguments (4 given)

In [134]: sparse.block_diag([DD2,DD2,DD2,DD2])
Out[134]: 
<16x16 sparse matrix of type '<type 'numpy.int32'>'
    with 40 stored elements in COOrdinate format>

这可能不是构建这种块对角数组的最快方法，但它是一个开始。

=================

查看sparse.block_mat 的代码，我推断确实如此：

In [145]: rows=[]
In [146]: for i in range(4):
    arow=[None]*4
    arow[i]=DD2
    rows.append(arow)
   .....:     

In [147]: rows
Out[147]: 
[[<4x4 sparse matrix of type '<type 'numpy.int32'>'
    with 10 stored elements (5 diagonals) in DIAgonal format>,
  None,
  None,
  None],
 [None,
  <4x4 sparse matrix of type '<type 'numpy.int32'>'
  ...
  None,
  <4x4 sparse matrix of type '<type 'numpy.int32'>'
    with 10 stored elements (5 diagonals) in DIAgonal format>]]

换句话说，rows 是None 的“矩阵”，DD2 沿对角线排列。然后它将这些传递给sparse.bmat。

In [148]: sparse.bmat(rows)
Out[148]: 
<16x16 sparse matrix of type '<type 'numpy.int32'>'
    with 40 stored elements in COOrdinate format>

bmat 反过来从所有输入矩阵的coo 格式中收集data,rows,cols，将它们加入主数组，并从中构建一个新的coo 矩阵。

所以另一种方法是直接构造这 3 个数组。

【讨论】：