块三对角矩阵蟒蛇

Question

我想从三个numpy.ndarray开始创建一个块三对角矩阵.在python中有没有(直接)方法呢？

先感谢您!

干杯

Answer 1

使用"常规"numpy数组,使用numpy.diag:

def tridiag(a, b, c, k1=-1, k2=0, k3=1):
    return np.diag(a, k1) + np.diag(b, k2) + np.diag(c, k3)

a = [1, 1]; b = [2, 2, 2]; c = [3, 3]
A = tridiag(a, b, c)

Answer 2

您还可以通过花式索引使用"常规"numpy数组执行此操作:

import numpy as np
data = np.zeros((10,10))
data[np.arange(5), np.arange(5)+2] = [5, 6, 7, 8, 9]
data[np.arange(3)+4, np.arange(3)] = [1, 2, 3]
print data

(你可以取代那些调用np.arange与np.r_如果你想成为更简洁.而不是如data[np.arange(3)+4, np.arange(3)]使用data[np.r_[:3]+4, np.r_[:3]])

这会产生:

[[0 0 5 0 0 0 0 0 0 0]
 [0 0 0 6 0 0 0 0 0 0]
 [0 0 0 0 7 0 0 0 0 0]
 [0 0 0 0 0 8 0 0 0 0]
 [1 0 0 0 0 0 9 0 0 0]
 [0 2 0 0 0 0 0 0 0 0]
 [0 0 3 0 0 0 0 0 0 0]
 [0 0 0 0 0 0 0 0 0 0]
 [0 0 0 0 0 0 0 0 0 0]
 [0 0 0 0 0 0 0 0 0 0]]

但是,如果您打算使用稀疏矩阵,请查看scipy.sparse.spdiags.(请注意,您必须在前面加上假数据到您的行值,如果你将数据放入一个对角线位置的正值(例如3点的位置4中的例子))

作为一个简单的例子:

import numpy as np
import scipy as sp
import scipy.sparse

diag_rows = np.array([[1, 1, 1, 1, 1, 1, 1],
                      [2, 2, 2, 2, 2, 2, 2],
                      [0, 0, 0, 0, 3, 3, 3]])
positions = [-3, 0, 4]
print sp.sparse.spdiags(diag_rows, positions, 10, 10).todense()

这会产生:

[[2 0 0 0 3 0 0 0 0 0]
 [0 2 0 0 0 3 0 0 0 0]
 [0 0 2 0 0 0 3 0 0 0]
 [1 0 0 2 0 0 0 0 0 0]
 [0 1 0 0 2 0 0 0 0 0]
 [0 0 1 0 0 2 0 0 0 0]
 [0 0 0 1 0 0 2 0 0 0]
 [0 0 0 0 1 0 0 0 0 0]
 [0 0 0 0 0 1 0 0 0 0]
 [0 0 0 0 0 0 1 0 0 0]]

Answer 3

@TheCorwoodRep的回答实际上可以在一行中完成。不需要单独的功能。

np.eye(3,3,k=-1) + np.eye(3,3)*2 + np.eye(3,3,k=1)*3

这产生：

array([[ 2.,  3.,  0.],
       [ 1.,  2.,  3.],
       [ 0.,  1.,  2.]])

Answer 4

使用该功能scipy.sparse.diags.

例:

from scipy.sparse import diags
import numpy as np
#
n = 10
k = np.array([np.ones(n-1),-2*np.ones(n),np.ones(n-1)])
offset = [-1,0,1]
A = diags(k,offset).toarray()

返回:

array([[-2.,  1.,  0.,  0.,  0.],
       [ 1., -2.,  1.,  0.,  0.],
       [ 0.,  1., -2.,  1.,  0.],
       [ 0.,  0.,  1., -2.,  1.],
       [ 0.,  0.,  0.,  1., -2.]])

Answer 5

由于三对角矩阵是一个稀疏矩阵，使用稀疏包可能是一个不错的选择，请参阅http://pysparse.sourceforge.net/spmatrix.html#matlab-implementation，其中有一些示例以及与 MATLAB 的比较，甚至...