2022年 12月 19日

python 对角阵_python-Numpy分区对角矩阵

我想像这样生成分区对角矩阵A

给出矩阵B

B = -np.diag(np.ones(n – 2), -1) – np.diag(np.ones(n – 2), 1) + 4 * np.diag(np.ones(n – 1))

例如,

有没有一种方法可以不使用循环?

抱歉,第一次错误地上传了矩阵A和B的图形.

解决方法:

您可以将构建块堆叠到查找表中,然后通过在其中建立索引来构建A:

>>> from scipy import sparse

>>>

>>> n = 5

>>> B = sparse.diags([-1, 4, -1], [-1, 0, 1], (n-1, n-1), dtype=int).A

>>> A = sparse.diags([1, 2, 1], [-1, 0, 1], (n-1, n-1), dtype=int).A

# 0 means 0 0 0 …,

# 1 means -I

# 2 means B

>>>

# next line builds the lookup table (using np.stack)

# does the lookup …[A]

# and flattens the resulting 4D array after swapping

# the middle axes; the swap reorders the entries from

# Vert, Horz, vert, horz to Vert, vert, Horz, horz

>>> A = np.stack([np.zeros_like(B), -np.identity(n-1, int), B])[A].swapaxes(1, 2).reshape((n-1)*(n-1), -1)

>>> A

array([[ 4, -1, 0, 0, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],

[-1, 4, -1, 0, 0, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0],

[ 0, -1, 4, -1, 0, 0, -1, 0, 0, 0, 0, 0, 0, 0, 0, 0],

[ 0, 0, -1, 4, 0, 0, 0, -1, 0, 0, 0, 0, 0, 0, 0, 0],

[-1, 0, 0, 0, 4, -1, 0, 0, -1, 0, 0, 0, 0, 0, 0, 0],

[ 0, -1, 0, 0, -1, 4, -1, 0, 0, -1, 0, 0, 0, 0, 0, 0],

[ 0, 0, -1, 0, 0, -1, 4, -1, 0, 0, -1, 0, 0, 0, 0, 0],

[ 0, 0, 0, -1, 0, 0, -1, 4, 0, 0, 0, -1, 0, 0, 0, 0],

[ 0, 0, 0, 0, -1, 0, 0, 0, 4, -1, 0, 0, -1, 0, 0, 0],

[ 0, 0, 0, 0, 0, -1, 0, 0, -1, 4, -1, 0, 0, -1, 0, 0],

[ 0, 0, 0, 0, 0, 0, -1, 0, 0, -1, 4, -1, 0, 0, -1, 0],

[ 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, -1, 4, 0, 0, 0, -1],

[ 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, 0, 4, -1, 0, 0],

[ 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, -1, 4, -1, 0],

[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, -1, 4, -1],

[ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, -1, 0, 0, -1, 4]])

请注意,稀疏构造函数仅出于其便利性而使用.稀疏矩阵将立即转换为密集矩阵(使用.A属性).

标签:python,numpy