spdiags（）如预期工作在Python

Question

在Matlab中/八度音功能不是，spdiags([-8037.500 50.000 -12.500], 0:2, 1, 51)给出以下输出：

(1, 1) -> -8037.5 
(1, 2) -> 50 
(1, 3) -> -12.500 

然而，当我使用在Python下面，它不会产生类似的结果作为在Matlab /八度：

​​

Python的spdiags()产生下面的输出，这是在缺少一号50和-12.5条款和第二指数：

array([[-8037.5,  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. ,  0. ,  0. ,  0. ,  0. , 
      0. ,  0. ,  0. ,  0. ,  0. ,  0. ,  0. , 
      0. ,  0. ]]) 

我看了一下this对一个类似问题的回答，但我不确定我哪里出错了。

编辑：

我想建立如下图所示的是由A_diag1，A_diag2和A_diag3矩阵A。按照答案中的建议，我定义了A_diag1和A_diag3。

import numpy as np 
import scipy as sp 
A_diag1 = np.tile(np.array([-8037.500, 50, -12.5]), (3,1)) 
A_diag2 = np.reshape(np.repeat([1250, -18505, 1250], 49), (3, 49)) 
A_diag3 = np.tile(np.array([12.5, -50, 8037.500]), (3,1)) 
A = np.concatenate((sp.sparse.spdiags(A_diag1, np.r_[0:2 + 1], 1, 51).toarray(), \ 
       sp.sparse.spdiags(A_diag2, np.r_[0:2 + 1], 49, 51).toarray(), \ 
       sp.sparse.spdiags(A_diag3, np.r_[48:50 + 1], 1, 51).toarray()), axis=0) 

然而，五显现细胞在过去的3行和A列显示为零/单数显示在下面的快照。我预计当前显示为零的突出显示的单元格不为零。 [您可以复制并粘贴上面的代码片段以再现A矩阵，从中可以获取如下所示的快照。]

EDIT2：

继使用sp.sparse.diags()按预期工作代码。与sp.sparse.spdiags不同，使用sp.sparse.diags()时，输入参数的结果形状（数组维数）必须在列表中。

import numpy as np 
import scipy as sp 
A_diag1 = np.array([[-8037.500], [50], [-12.5]]) 
A_diag2 = np.reshape(np.repeat([1250, -18505, 1250], 49), (3, 49)) 
A_diag3 = np.array([[12.5], [-50], [8037.500]]) 
A = np.concatenate((sp.sparse.diags(A_diag1, np.arange(0, 2 + 1), [1, 51]).toarray(), \ 
sp.sparse.diags(A_diag2, np.arange(0, 2 + 1), [49, 51]).toarray(), \ 
sp.sparse.diags(A_diag3, np.arange(48, 50 + 1), [1, 51]).toarray()), axis=0) 

Answer 1

这使得一个稀疏矩阵（51.1），具有值向下每一行：

In [5]: sparse.spdiags(data,[0,-1,-2],51,1) 
Out[5]: 
<51x1 sparse matrix of type '<class 'numpy.float64'>' 
    with 3 stored elements (3 diagonals) in DIAgonal format> 
In [6]: print(_) 
    (0, 0) -8037.5 
    (1, 0) 50.0 
    (2, 0) -12.5 

注意，spdiags定义：

数据：存储array_like 矩阵对角线行向明智的

稀疏diagonal format将其数据存储在矩阵中，其中的一部分可以是“屏幕外”。所以使用起来有点棘手。我通常使用coo类型的输入创建矩阵。

In [27]: M =sparse.spdiags(data,[0,-1,-2],3,3) 
In [28]: M.A 
Out[28]: 
array([[-8037.5,  0. ,  0. ], 
     [ 50. ,  0. ,  0. ], 
     [ -12.5,  0. ,  0. ]]) 
In [29]: M.data 
Out[29]: 
array([[-8037.5], 
     [ 50. ], 
     [ -12.5]]) 
In [30]: M.offsets 
Out[30]: array([ 0, -1, -2], dtype=int32) 

你想要的是它的转置（也许）

In [32]: Mt = M.T 
In [33]: Mt.A 
Out[33]: 
array([[-8037.5, 50. , -12.5], 
     [ 0. ,  0. ,  0. ], 
     [ 0. ,  0. ,  0. ]]) 
In [34]: Mt.data 
Out[34]: 
array([[-8037.5,  0. ,  0. ], 
     [ 0. , 50. ,  0. ], 
     [ 0. ,  0. , -12.5]]) 
In [35]: Mt.offsets 
Out[35]: array([0, 1, 2], dtype=int32) 

因此，我们可以重新Mt有：

sparse.spdiags(Mt.data, Mt.offsets, 3,3) 

---

如果我救了八度矩阵并加载它，我得到：

In [40]: loadmat('diags') 
Out[40]: 
{'__globals__': [], 
'__header__': b'MATLAB 5.0 MAT-file, written by Octave 4.0.0, 2017-10-19 01:24:58 UTC', 
'__version__': '1.0', 
'x': <1x51 sparse matrix of type '<class 'numpy.float64'>' 
    with 3 stored elements in Compressed Sparse Column format>} 
In [42]: X=_['x'] 
In [43]: print(X) 
    (0, 0) -8037.5 
    (0, 1) 50.0 
    (0, 2) -12.5 

如果我把它转换为dia格式我得到类似Mt：

In [48]: sparse.dia_matrix(X) 
Out[48]: 
<1x51 sparse matrix of type '<class 'numpy.float64'>' 
    with 3 stored elements (3 diagonals) in DIAgonal format> 
In [49]: print(_) 
    (0, 0) -8037.5 
    (0, 1) 50.0 
    (0, 2) -12.5 
In [50]: _.data, _.offsets 
Out[50]: 
(array([[-8037.5,  0. ,  0. ], 
     [ 0. , 50. ,  0. ], 
     [ 0. ,  0. , -12.5]]), array([0, 1, 2])) 

的sparse.diags功能可能更直观：

In [92]: sparse.diags(data, [0,1,2],(1,3)) 
Out[92]: 
<1x3 sparse matrix of type '<class 'numpy.float64'>' 
    with 3 stored elements (3 diagonals) in DIAgonal format> 
In [93]: _.A 
Out[93]: array([[-8037.5, 50. , -12.5]]) 
In [94]: print(__) 
    (0, 0) -8037.5 
    (0, 1) 50.0 
    (0, 2) -12.5 

---

In [56]: sp1 = sparse.spdiags(A_diag1, np.r_[0:2 + 1], 1, 51) 
In [57]: sp2 = sparse.spdiags(A_diag2, np.r_[0:2 + 1], 49, 51) 
In [58]: sp3 = sparse.spdiags(A_diag3, np.r_[48:50 + 1], 1, 51) 

（该r_表达式也可以是np.arange(0,3)和np.arange(48,51)）

这些可以sparse.vstack加盟（它结合了coo格式属性）

In [69]: B = sparse.vstack((sp1,sp2,sp3)) 
    In [72]: B 
    Out[72]: 
    <51x51 sparse matrix of type '<class 'numpy.float64'>' 
     with 147 stored elements in COOrdinate format> 

In [75]: B.tocsr()[45:, 46:].A 
Out[75]: 
array([[ 1250.,  0.,  0.,  0.,  0.], 
     [-18505., 1250.,  0.,  0.,  0.], 
     [ 1250., -18505., 1250.,  0.,  0.], 
     [  0., 1250., -18505.,  0.,  0.], 
     [  0.,  0., 1250.,  0.,  0.], 
     [  0.,  0.,  0.,  0.,  0.]]) 

快照匹配。 （我仍然需要弄清楚你正在创建什么）。

sparse.spdiags(data, diags, m, n)是调用sparse.dia_matrix((data, diags), shape=(m,n))

再回到sparse.diags，如果你想3个对角线，每个充满了来自data一个值，我们可以用另一种方式：

​​

所以sp1将有看起来像

In [117]: B.tocsr()[0,:].todia().data 
Out[117]: 
array([[-8037.5,  0. ,  0. ], 
     [ 0. , 50. ,  0. ], 
     [ 0. ,  0. , -12.5]]) 

Answer 2

我对你的观察没有解释（没有太大的MATLAB用户的，但我可以证实，倍频是做什么像你说的），但以下SciPy的的example-usage，则可以使用达到这一结果：

import numpy as np 
import scipy.sparse as sp 

data = np.tile(np.array([-8037.5, 50., -12.5]), (3,1)) 
x = sp.spdiags(data, np.arange(3), 1, 51) 

print(x) 

输出：

(0, 0)  -8037.5 
(0, 1)  50.0 
(0, 2)  -12.5 

瓦片步骤构建：

[[-8037.5 50. -12.5] 
[-8037.5 50. -12.5] 
[-8037.5 50. -12.5]] 

当然，一切都是基于0索引的。