Python: 创建(稀疏)堆叠对角块矩阵

huangapple go评论46阅读模式
英文:

Python: create (sparse) stacked diagonal block matrix

问题

我需要创建一个具有以下形式的矩阵:

M=[
[a1, 0, 0],
[0, b1, 0],
[0, 0, c1],
[a2, 0, 0],
[0, b2, 0],
[0, 0, c2],
[a3, 0, 0],
[0, b3, 0],
[0, 0, c3],
...]

其中a(i),b(i)和c(i)是[1xp]块。生成的矩阵M具有[3m x 3p]的形式。我已经以3个矩阵[m x p]的形式获得了输入数据:

A = [[a1.T, a2.T, a3.T, ...]].T
B = [[b1.T, b2.T, b3.T, ...]].T
C = [[c1.T, c2.T, c3.T, ...]].T

我应该如何创建矩阵M?理想情况下,它会使用scipy.sparse库创建稀疏矩阵,但即使使用numpy创建稠密矩阵也让我很困扰。在这种情况下是否无法避免使用循环或至少使用列表推导式?

英文:

I need to create a matrix with the form

M=[
[a1, 0, 0],
[0, b1, 0],
[0, 0, c1],
[a2, 0, 0],
[0, b2, 0],
[0, 0, c2],
[a3, 0, 0],
[0, b3, 0],
[0, 0, c3],
...]

where a(i), b(i) and c(i) are [1xp] blocks. The resulting matrix M has the form [3m x 3p]. I am given the input data in the form of 3 matrices [m x p]:

A = [[a1.T, a2.T, a3.T, ...]].T
B = [[b1.T, b2.T, b3.T, ...]].T
C = [[c1.T, c2.T, c3.T, ...]].T

How can I create the matrix M? Ideally it would be sparse using the scipy.sparse library but I am even struggling creating it as a dense matrix using numpy. Is there no way around a loop or at least list comprehension in this case?

答案1

得分: 1

No need to make it complicated. For your scale, the following executes in less than a second.

import numpy as np
import scipy.sparse
from numpy.random import default_rng

rand = default_rng(seed=0)
m = 70,000
p = 20
abc = rand.random((3, m, p))
M_dense = np.zeros((m, 3, 3*p))
for i in range(3):
    M_dense[:, i, i*p:(i+1)*p] = abc[i, ...]
M_sparse = scipy.sparse.csr_matrix(M_dense.reshape((-1, 3*p)))
print(M_sparse.shape)
(210,000, 60)

Far better, though, is to construct the sparse matrix directly. Note the permuted shape of abc.

abc = rand.random((m, 3, p))
data = abc.ravel()
indices = np.tile(np.arange(3*p), m)
indptr = np.arange(0, data.size+1, p)
M_sparse = scipy.sparse.csr_matrix((data, indices, indptr))
英文:

No need to make it complicated. For your scale, the following executes in less than a second.

import numpy as np
import scipy.sparse
from numpy.random import default_rng

rand = default_rng(seed=0)
m = 70_000
p = 20
abc = rand.random((3, m, p))
M_dense = np.zeros((m, 3, 3*p))
for i in range(3):
    M_dense[:, i, i*p:(i+1)*p] = abc[i, ...]
M_sparse = scipy.sparse.csr_matrix(M_dense.reshape((-1, 3*p)))
print(M_sparse.shape)
(210000, 60)

Far better, though, is to construct the sparse matrix directly. Note the permuted shape of abc.

abc = rand.random((m, 3, p))
data = abc.ravel()
indices = np.tile(np.arange(3*p), m)
indptr = np.arange(0, data.size+1, p)
M_sparse = scipy.sparse.csr_matrix((data, indices, indptr))

huangapple
  • 本文由 发表于 2023年2月19日 01:56:54
  • 转载请务必保留本文链接:https://go.coder-hub.com/75495303.html
匿名

发表评论

匿名网友

:?: :razz: :sad: :evil: :!: :smile: :oops: :grin: :eek: :shock: :???: :cool: :lol: :mad: :twisted: :roll: :wink: :idea: :arrow: :neutral: :cry: :mrgreen:

确定