1 2 #include <petscmat.h> 3 #include <petsc/private/petscimpl.h> 4 5 /* 6 MatDumpSPAI - Dumps a PETSc matrix to a file in an ASCII format 7 suitable for the SPAI code of Stephen Barnard to solve. This routine 8 is simply here to allow testing of matrices directly with the SPAI 9 code, rather then through the PETSc interface. 10 11 */ 12 PetscErrorCode MatDumpSPAI(Mat A, FILE *file) 13 { 14 PetscMPIInt size; 15 PetscInt n; 16 MPI_Comm comm; 17 18 PetscFunctionBegin; 19 PetscValidHeaderSpecific(A, MAT_CLASSID, 1); 20 PetscAssertPointer(file, 2); 21 PetscCall(PetscObjectGetComm((PetscObject)A, &comm)); 22 PetscCallMPI(MPI_Comm_size(comm, &size)); 23 PetscCheck(size == 1, comm, PETSC_ERR_SUP, "Only single processor dumps"); 24 PetscCall(MatGetSize(A, &n, &n)); 25 /* print the matrix */ 26 fprintf(file, "%" PetscInt_FMT "\n", n); 27 for (PetscInt i = 0; i < n; i++) { 28 const PetscInt *cols; 29 const PetscScalar *vals; 30 PetscInt nz; 31 32 PetscCall(MatGetRow(A, i, &nz, &cols, &vals)); 33 for (PetscInt j = 0; j < nz; j++) fprintf(file, "%" PetscInt_FMT " %d" PetscInt_FMT " %16.14e\n", i + 1, cols[j] + 1, vals[j]); 34 PetscCall(MatRestoreRow(A, i, &nz, &cols, &vals)); 35 } 36 PetscFunctionReturn(PETSC_SUCCESS); 37 } 38 39 PetscErrorCode VecDumpSPAI(Vec b, FILE *file) 40 { 41 PetscInt n; 42 const PetscScalar *array; 43 44 PetscFunctionBegin; 45 PetscValidHeaderSpecific(b, VEC_CLASSID, 1); 46 PetscAssertPointer(file, 2); 47 PetscCall(VecGetSize(b, &n)); 48 PetscCall(VecGetArrayRead(b, &array)); 49 fprintf(file, "%" PetscInt_FMT "\n", n); 50 for (PetscInt i = 0; i < n; i++) fprintf(file, "%" PetscInt_FMT " %16.14e\n", i + 1, array[i]); 51 PetscCall(VecRestoreArrayRead(b, &array)); 52 PetscFunctionReturn(PETSC_SUCCESS); 53 } 54