scipy dia_matrix 的 Python 示例

Question

运行这段代码

d0  = np.ones(N)
dp1 = np.ones(N - 1)
dm1 = np.ones(N - 1)

diag = [[d0],[dp1],[dm1]]
offsets = [0,1,-1]

A = dia_matrix( (diag,offsets), shape=(N,N), dtype=float)

我收到以下错误:

File "/usr/local/lib/python2.7/dist-packages/scipy/sparse/dia.py", line 109, in __init__
self.data = np.atleast_2d(np.array(arg1[0], dtype=dtype, copy=copy))

ValueError:使用序列设置数组元素。

我不明白我做错了什么!有人可以为我提供一个正确的例子来做我想做的事吗？

Answer 1

当 dia_matrix 的第一个参数时具有 (data, offsets) 的形式，data 应该是一个二维数组，data 的每一行都有一条对角线矩阵。因为 data 是一个矩形矩阵，所以 data 中的一些元素被忽略了。次对角线是“左对齐”，上对角线是“右对齐”。 (具体来说，data 和稀疏矩阵 A 之间的映射是 data[i,j] == A[j - offsets[i], j].) 例如，考虑以下将用于创建 5x5 矩阵的内容:

In [28]: data
Out[28]: 
array([[ 1,  2,  3,  4,  5],
       [ 6,  7,  8,  9, 10],
       [11, 12, 13, 14, 15]])

In [29]: offsets
Out[29]: [0, 1, -1]

data 包含三个对角线。由于 offset[0] 为 0，因此 data 的第 0 行包含主对角线。此行中的所有 5 个元素都用于矩阵中。 offset[1] 为 1，因此 data[1] 中的数据成为第一个超对角线。仅使用值 [7, 8, 9, 10]；第一个值 6 被忽略。类似地，data 的第三行给出第一个次对角线，并且仅使用值 [11, 12, 13, 14]。

In [30]: a = dia_matrix((data, offsets), shape=(5, 5))

In [31]: a.A
Out[31]: 
array([[ 1,  7,  0,  0,  0],
       [11,  2,  8,  0,  0],
       [ 0, 12,  3,  9,  0],
       [ 0,  0, 13,  4, 10],
       [ 0,  0,  0, 14,  5]])

您的示例可以重写如下:

In [32]: N = 5

In [33]: data = np.ones((3, 5))

In [34]: A = dia_matrix((data, offsets), shape=(N, N), dtype=float)

In [35]: A.A
Out[35]: 
array([[ 1.,  1.,  0.,  0.,  0.],
       [ 1.,  1.,  1.,  0.,  0.],
       [ 0.,  1.,  1.,  1.,  0.],
       [ 0.,  0.,  1.,  1.,  1.],
       [ 0.,  0.,  0.,  1.,  1.]])

dia_matrix 文档字符串还有另一个示例。

或者，您可以使用 scipy.sparse.diags创建矩阵。如果您已经有了生成“正确”大小的对角线的代码，这将很有用。使用 diags，您不必创建矩形 data 矩阵。例如，

In [104]: from scipy.sparse import diags

In [105]: d0 = ones(n)

In [106]: dp1 = np.ones(N - 1)

In [107]: dm1 = np.ones(N - 1)

In [108]: d = [d0, dp1, dm1]

In [109]: B = diags(d, offsets, dtype=float)

In [110]: B.A
Out[110]: 
array([[ 1.,  1.,  0.,  0.,  0.],
       [ 1.,  1.,  1.,  0.,  0.],
       [ 0.,  1.,  1.,  1.,  0.],
       [ 0.,  0.,  1.,  1.,  1.],
       [ 0.,  0.,  0.,  1.,  1.]])