【问题标题】:Read sparse matrix from ascii file into python将 ascii 文件中的稀疏矩阵读入 python
【发布时间】:2019-08-09 15:53:50
【问题描述】:

我有一个使用并行 PETSc 稀疏矩阵格式 mpiaij 的 Fortran 代码。

我想对这些矩阵做一些分析,所以我想将它们读入 python。

我尝试了 Fortran 中的二进制输出和 petsc4py 中的二进制输入,但显然它们不兼容。 Petsc HDF5 输出会创建不可读的 HDF5 文件,所以我现在只能使用 ASCII 格式。

在 ascii 中,矩阵如下所示:

Mat Object: 48 MPI processes
  type: mpiaij
row 0: (0, 0.934865)  (1, 0.00582401)  (2, -0.00125881)  (3, 0.000157352)  (10, 0.0212704)  (11, -9.37151e-05)  (12, 7.77296e-06)  (13, 1.15276e-06)  (20, -0.00457321)  (21, 9.31045e-06)  (22, -1.37541e-07)  (23, -3.00994e-07)  (30, 0.000571716)  (31, 5.82622e-07)  (32, -2.27908e-07)  (33, 4.55904e-08)  (3410, 0.0005718)  (3411, 3.14914e-06)  (3412, -5.83246e-07)  (3413, 5.58045e-08)  (3420, -0.00457491)  (3421, -3.91645e-05)  (3422, 6.62677e-06)  (3423, -5.10165e-07)  (3430, 0.0212818)  (3431, 0.000230778)  (3432, -3.75686e-05)  (3433, 2.57173e-06) 
row 1: (...)

有没有一种优雅的方法可以将它解析成 python?

【问题讨论】:

    标签: python sparse-matrix petsc


    【解决方案1】:

    我不熟悉 PETSc 或其矩阵格式,但鉴于示例 ASCII 格式,当然可以将其转换为 Python 中的任何其他矩阵格式。我假设该文件包含每个非零行的一行,并且每行中的数字对是列索引和相应的数字。对吗?

    您认为“一种优雅的方式”是个人意见,对于 Stack Overflow 来说并不是一个真正有效的问题,但我可以尝试为您指出可行的解决方案的正确方向。

    首先,在不了解所有细节的情况下,在我看来,正确的问题应该是“为什么 Fortran 中的二进制输出和 petsc4py 中的二进制输入不兼容?”如果你能解决这个问题,那可能是最好的解决方案。如果我没记错的话,Fortran 代码支持不同的字节顺序,默认情况下可能使用大端格式,而 Python 通常使用小端格式。也许您可以在其中一个库函数中指定字节顺序,或者您可以在必要时手动转换字节顺序。这可能是您首先要研究的内容。

    作为一种解决方法,您可以在 Python 中解析 A​​SCII 格式以进行进一步处理。我假设您已经搜索了现有库并且找不到任何库,因此您需要编写一些自定义代码。根据您的需要,“不错”的解决方案将使用正则表达式,但一种快速而简单的方法是使用标准字符串方法和 eval() 函数,因为 ASCII 格式已经非常类似于 Python 语法:-)

    注意:如果您信任输入文件,请仅使用eval() 函数,因为它容易受到代码注入攻击!对于个人使用,这通常不是问题。

    我在下面提供了一些示例代码。这完成了基本的输入处理。你想对数据做什么取决于你,所以你需要自己完成代码。此示例代码仅打印数字。

    def read_mpiaij(file):
        lines = file.read().splitlines()
        assert 'Mat Object: ' in lines[0]
        assert lines[1] == '  type: mpiaij'
        for line in lines[2:]:
            parts = line.split(': ')
            assert len(parts) == 2
            assert parts[0].startswith('row ')
    
            row_index = int(parts[0][4:])
            row_contents = eval(parts[1].replace(')  (', '), ('))
    
            # Here you have the row_index and a tuple of (column_index, value)
            # pairs that specify the non-zero contents. You could process this
            # depending on your needs, e.g. store the values in an array.
            for (col_index, value) in row_contents:
                print('row %d, col %d: %s' % (row_index, col_index, value))
                # TODO: Implement real code here.
                # You probably want to do something like:
                # data[row_index][col_index] = value
    
    
    def main():
        with open('input.txt', 'rt', encoding='ascii') as file:
            read_mpiaij(file)
    
    
    if __name__ == '__main__':
        main()
    

    输出:

    row 0, col 0: 0.934865
    row 0, col 1: 0.00582401
    row 0, col 2: -0.00125881
    row 0, col 3: 0.000157352
    row 0, col 10: 0.0212704
    row 0, col 11: -9.37151e-05
    row 0, col 12: 7.77296e-06
    row 0, col 13: 1.15276e-06
    row 0, col 20: -0.00457321
    row 0, col 21: 9.31045e-06
    row 0, col 22: -1.37541e-07
    row 0, col 23: -3.00994e-07
    row 0, col 30: 0.000571716
    row 0, col 31: 5.82622e-07
    row 0, col 32: -2.27908e-07
    row 0, col 33: 4.55904e-08
    row 0, col 3410: 0.0005718
    row 0, col 3411: 3.14914e-06
    row 0, col 3412: -5.83246e-07
    row 0, col 3413: 5.58045e-08
    row 0, col 3420: -0.00457491
    row 0, col 3421: -3.91645e-05
    row 0, col 3422: 6.62677e-06
    row 0, col 3423: -5.10165e-07
    row 0, col 3430: 0.0212818
    row 0, col 3431: 0.000230778
    row 0, col 3432: -3.75686e-05
    row 0, col 3433: 2.57173e-06
    ...
    

    【讨论】:

      【解决方案2】:

      正则表达式是你的朋友。像这样的东西怎么样:

      for recnum, rec in enumerate(fh.readlines()):
          mat = re.match(r'row\s*(\d+):\s*(.*)', rec)
          if (not mat): raise IOError("Bad data at rec %d." % (recnum))
          rowNum = int(mat.group(1))
          rest = mat.group(2)
          lastColNum = -1
          for col in re.finditer(r'\(\d+),\s*(\d+\.\d*\)', rest):
              colNum = int(mat.group(1))
              if (colNum <= lastColNum):
                  raise KeyError("colNum out of order at rec %d." % (colNum, recNum))
              value = float(mat.group(2))
              # save cell, like via numpy tbl[rowNum, colNum] = value
      

      我假设每行中的列项是按顺序排列的。如果没有,或者有其他约束(例如,如果值必须在 0.0...1.0 中,这在您的示例中似乎是正确的),您当然可以进行调整。值得检查数据,因为数据很少像人们希望的那样干净......

      【讨论】:

        猜你喜欢
        • 2013-11-05
        • 2014-11-01
        • 2023-04-10
        • 2018-01-19
        • 2015-04-26
        • 2023-03-17
        • 1970-01-01
        • 2020-11-27
        • 2011-11-20
        相关资源
        最近更新 更多