I am trying to read a sparse PETSc matrix saved from a Fortran code into a binary file like so:
CALL PetscViewerBinaryOpen(PETSC_COMM_SELF, TRIM(ADJUSTL(filename)), FILE_MODE_WRITE, viewer2, ier)
CALL MatView(aa_symmetric, viewer2, ier)
CALL PetscViewerDestroy(viewer2, ier)
When I run the Fortran code single-threaded, all works fine and I can read the matrix without any issues using the following Python code:
import petsc4py
from petsc4py import PETSc
# Initialize PETSc
petsc4py.init(sys.argv)
# Create a viewer for reading the binary file
viewer = PETSc.Viewer().createBinary(filename, mode='r', comm=PETSc.COMM_WORLD)
# Create a matrix and load data from the binary file
A_petsc = PETSc.Mat().create(comm=PETSc.COMM_WORLD)
A_petsc.setType(PETSc.Mat.Type.AIJ)
A_petsc.setFromOptions()
A_petsc.load(viewer)
# Finalize PETSc
PETSc.Finalize()
However, when I run the Fortran code on more processors ("mpirun -n 2 my_exe"), I get the following error on the Python side:
Traceback (most recent call last):
File "/home/Python/test_matrixImport_binary.py", line 80, in <module>
A_petsc.load(viewer) File "petsc4py/PETSc/Mat.pyx", line 2025, in
petsc4py.PETSc.Mat.load petsc4py.PETSc.Error: error code 79 [0] MatLoad() at
/home/lib/petsc-3.21.0/src/mat/interface/matrix.c:1344 [0] MatLoad_SeqAIJ() at
/home/lib/petsc-3.21.0/src/mat/impls/aij/seq/aij.c:5091 [0] MatLoad_SeqAIJ_Binary() at
/home/lib/petsc-3.21.0/src/mat/impls/aij/seq/aij.c:5142 [0] Unexpected data in file [0]
Inconsistent matrix data in file: nonzeros = 460, sum-row-lengths = 761
How can I fix this?
It seems that the issue was on the export side, not on the PETSc side. I originally wrote the export routine for symmetric matrices but then applied it to one that's not symmetric, leading to confused addressing in the binary file. I hadn't caught it earlier because I was using my own ASCII import routine on the Python side that didn't care about symmetry. The snippet in the question should work fine. I am now using the following two functions that I think are a bit cleaner:
import numpy as np
from scipy.sparse import coo_matrix, csr_matrix
from petsc4py import PETSc
def readMat(filename, toScipy=False):
viewer = PETSc.Viewer().createBinary(filename, mode='r', comm=PETSc.COMM_WORLD)
A_petsc = PETSc.Mat().create(comm=PETSc.COMM_WORLD)
A_petsc.setType(PETSc.Mat.Type.AIJ)
A_petsc.setFromOptions()
A_petsc.load(viewer)
if toScipy:
I, J, V = A_petsc.getValuesCSR()
return csr_matrix((V, J, I)).tocoo()
else:
return A_petsc
def readVec(filename, toNumpy=False):
viewer = PETSc.Viewer().createBinary(filename, mode='r', comm=PETSc.COMM_WORLD)
v_petsc = PETSc.Vec().create(comm=PETSc.COMM_WORLD)
v_petsc.setFromOptions()
v_petsc.load(viewer)
if toNumpy:
return np.array(v_petsc)
else:
return v_petsc