Matlab 稀疏函数与Python Scipy.sparse 库答案

【问题标题】：Matlab Sparse function v.s. Python Scipy.sparse libraryMatlab 稀疏函数与Python Scipy.sparse 库
【发布时间】：2020-05-11 22:03:51
【问题描述】：

我正在将 Matlab 代码翻译成 Python 环境，现在在 Matlab 中处理稀疏函数。我知道我们有一个名为 scipy.sparse 的库，其中一个名为 csc_matrix((data, (row_ind, col_ind)), [shape=(M, N)])。

但是，当我检查从 scipy.sparse 库计算的那些时，它们与 Matlab 中的不匹配。我有大小为 1693872 的大数据，称为 Ig (1693872,)、Jg(1693872,) 和 K_dummy (1693872,)，其中 K_dummy(Ig(i),Jg(i)) = K_dummy(i)。

我已经用 Matlab 检查了所有变量 Ig、Jg、K_dummy 并且完全匹配。你们知道我必须考虑其他方面吗？

这是我在 python 和 Matlab 中的示例代码，分别作为参考：

K = csc_matrix((K_dummy.flatten('F'),(Ig.flatten('F')-1,Jg.flatten('F')-1)),shape=(noDofs,noDofs))

K = sparse(Ig(:),Jg(:),K_dummy_python(:),noDofs,noDofs);

其中K_dummy 是 (18, 18, 5228) 数组，Ig 是 (324, 5228) 数组，Jg 是 (324, 5228) 数组，noDofs 是一个 int 变量，为 42442。

【问题讨论】：

如果示例很小且我们可以在自己的机器上重现值，则最容易提供帮助。虽然几年前我曾在 MATLAB 上工作过，但我还没有尝试过并排创作（使用 Octave）。我已经通过 .mat 文件导入 Octave 创建稀疏。
我很想在 Octave 中创建一个样本，然后将密集和稀疏的值传输到 scipy（通过 io.loadmat），但是制作一些能够捕捉形状和顺序复杂性的东西正在推进比我现在要做的更多的工作。有什么区别？稀疏性（非零位置）是否正确？只是非零值的差异？
它们之间的区别在于存储在稀疏矩阵中的非零元素的数量......而且，我检查了所有扁平化的索引，Ig、Jg 和所有正确的。
应该有 1693872 个非零元素（M.nnz、nnz(M)）。可能的例外 - K_dummy 有一些 0 值，其中一个正在删除，另一个没有。或者有一些重复，一个是求和，另一个不是。 scipy 有一个 M.eliminate_zeros() 方法可以在处理留下零时使用。对重复项求和是两者的常态。

标签： python matlab scipy sparse-matrix

【解决方案1】：

八度：

>> data=[1,2;3,4];
>> I=[1,2,3,4];
>> J=[2,3,4,2];
>> M = sparse(I(:),J(:),data(:))
M =

Compressed Column Sparse (rows = 4, cols = 4, nnz = 4 [25%])

  (1, 2) ->  1
  (4, 2) ->  4
  (2, 3) ->  3
  (3, 4) ->  2

>> full(M)
ans =

   0   1   0   0
   0   0   3   0
   0   0   0   2
   0   4   0   0

>> save -7 sparse1.mat data I J M

在numpy 与scipy.io 和scipy.sparse：

In [57]: d = loadmat('sparse1.mat')                                                              
In [58]: d                                                                                       
Out[58]: 
{'__header__': b'MATLAB 5.0 MAT-file, written by Octave 4.2.2, 2020-01-26 19:36:04 UTC',
 '__version__': '1.0',
 '__globals__': [],
 'data': array([[1., 2.],
        [3., 4.]]),
 'I': array([[1., 2., 3., 4.]]),
 'J': array([[2., 3., 4., 2.]]),
 'M': <4x4 sparse matrix of type '<class 'numpy.float64'>'
    with 4 stored elements in Compressed Sparse Column format>}

并查看稀疏矩阵：

In [59]: d['M'].A                                                                                
Out[59]: 
array([[0., 1., 0., 0.],
       [0., 0., 3., 0.],
       [0., 0., 0., 2.],
       [0., 4., 0., 0.]])
In [60]: print(d['M'])                                                                           
  (0, 1)    1.0
  (3, 1)    4.0
  (1, 2)    3.0
  (2, 3)    2.0

然后重新创建矩阵

In [61]: M1 = sparse.csc_matrix((d['data'].flatten('F'), 
                (d['I'].astype(int).flatten('F')-1, 
                 d['J'].astype(int).flatten('F')-1)))                                                        
In [62]: M1                                                                                      
Out[62]: 
<4x4 sparse matrix of type '<class 'numpy.float64'>'
    with 4 stored elements in Compressed Sparse Column format>
In [63]: M1.A                                                                                    
Out[63]: 
array([[0., 1., 0., 0.],
       [0., 0., 3., 0.],
       [0., 0., 0., 2.],
       [0., 4., 0., 0.]])
In [64]: print(M1)                                                                               
  (0, 1)    1.0
  (3, 1)    4.0
  (1, 2)    3.0
  (2, 3)    2.0

有了更新的numpy/scipy，我不得不将索引转换为整数（他们对浮点数作为索引变得更加挑剔了）。但是扁平化似乎效果很好。

你在所有这些变量中都有一个额外的维度，我懒得重新创建。正确展平可能会导致您的问题。

【讨论】：