这是我的第一个SO问题。如果我可以问得更好，请告诉我：）
我正试图找到一种方法，将稀疏矩阵列表拼接成一个更大的块矩阵。
我有一个python代码，它生成正方形稀疏矩阵的列表，一个矩阵接一个矩阵。在伪代码中：
Lx = [Lx1, Lx1, ... Lxn]
Ly = [Ly1, Ly2, ... Lyn]
Lz = [Lz1, Lz2, ... Lzn]
由于每个单独的Lx1、Lx2等矩阵都是按顺序计算的，所以它们被附加到一个列表中——我找不到一种方法来"动态"填充一个类似数组的对象。
我正在优化速度，瓶颈的特点是逐项计算笛卡尔乘积，类似于伪代码：
M += J[i,j] * [ Lxi *Lxj + Lyi*Lyj + Lzi*Lzj ]
对于0<=的所有组合i、 j<=n.（J是一个n维的数字平方矩阵）。
似乎通过（伪代码）一步计算所有笛卡尔乘积来对其进行矢量化：
L = [ [Lx1, Lx2, ...Lxn],
[Ly1, Ly2, ...Lyn],
[Lz1, Lz2, ...Lzn] ]
product = L.T * L
会更快。然而，像np.bmp、np.vstack、np.hstack这样的选项似乎需要数组作为输入，而我有列表。
有没有一种方法可以有效地将三个矩阵列表拼接成一个块？或者，有没有一种方法可以生成一个稀疏矩阵阵列，一次生成一个元素，然后将它们np.vstack在一起？
参考：类似的MATLAB代码，用于计算n自旋NMR模拟的哈密顿矩阵，可以在这里找到：
http://spindynamics.org/Spin-Dynamics---Part-II---Lecture-06.php
这是scipy.sparse.bmat:
L = scipy.sparse.bmat([Lx, Ly, Lz], format='csc')
LT = scipy.sparse.bmat(zip(Lx, Ly, Lz), format='csr') # Not equivalent to L.T
product = LT * L
我有一个"矢量化"的解决方案，但它的速度几乎是原始代码的两倍。根据kernprof测试，上面显示的瓶颈和下面最后一行显示的最终点积都需要大约95%的计算时间。
# Create the matrix of column vectors from these lists
L_column = bmat([Lx, Ly, Lz], format='csc')
# Create the matrix of row vectors (via a transpose of matrix with
# transposed blocks)
Lx_trans = [x.T for x in Lx]
Ly_trans = [y.T for y in Ly]
Lz_trans = [z.T for z in Lz]
L_row = bmat([Lx_trans, Ly_trans, Lz_trans], format='csr').T
product = L_row * L_column
使用稀疏矩阵和数组，我能够通过而不是将速度提高十倍。
Lx = np.empty((1, nspins), dtype='object')
Ly = np.empty((1, nspins), dtype='object')
Lz = np.empty((1, nspins), dtype='object')
这些数组在生成时由单独的Lx数组（以前是稀疏矩阵）填充。使用阵列结构允许转置和笛卡尔乘积按需执行：
Lcol = np.vstack((Lx, Ly, Lz)).real
Lrow = Lcol.T # As opposed to sparse version of code, this works!
Lproduct = np.dot(Lrow, Lcol)
单独的Lx[n]矩阵仍然是"捆绑的"，因此乘积是一个nxn矩阵。这意味着n x n J阵列与L乘积的原位乘法工作：
scalars = np.multiply(J, Lproduct)
然后将每个矩阵元素添加到最终的哈密尔顿矩阵：
for n in range(nspins):
for m in range(nspins):
M += scalars[n, k].real