1: #include <petscmat.h>
  2: #include <petsc/private/petscimpl.h>

  4: /*
  5:      MatDumpSPAI - Dumps a PETSc matrix to a file in an ASCII format
  6:   suitable for the SPAI code of Stephen Barnard to solve. This routine
  7:   is simply here to allow testing of matrices directly with the SPAI
  8:   code, rather than through the PETSc interface.

 10: */
 11: PetscErrorCode MatDumpSPAI(Mat A, FILE *file)
 12: {
 13:   PetscMPIInt size;
 14:   PetscInt    n;
 15:   MPI_Comm    comm;

 17:   PetscFunctionBegin;
 19:   PetscAssertPointer(file, 2);
 20:   PetscCall(PetscObjectGetComm((PetscObject)A, &comm));
 21:   PetscCallMPI(MPI_Comm_size(comm, &size));
 22:   PetscCheck(size == 1, comm, PETSC_ERR_SUP, "Only single processor dumps");
 23:   PetscCall(MatGetSize(A, &n, &n));
 24:   /* print the matrix */
 25:   fprintf(file, "%" PetscInt_FMT "\n", n);
 26:   for (PetscInt i = 0; i < n; i++) {
 27:     const PetscInt    *cols;
 28:     const PetscScalar *vals;
 29:     PetscInt           nz;

 31:     PetscCall(MatGetRow(A, i, &nz, &cols, &vals));
 32:     for (PetscInt j = 0; j < nz; j++) fprintf(file, "%" PetscInt_FMT " %d" PetscInt_FMT " %16.14e\n", i + 1, cols[j] + 1, vals[j]);
 33:     PetscCall(MatRestoreRow(A, i, &nz, &cols, &vals));
 34:   }
 35:   PetscFunctionReturn(PETSC_SUCCESS);
 36: }

 38: PetscErrorCode VecDumpSPAI(Vec b, FILE *file)
 39: {
 40:   PetscInt           n;
 41:   const PetscScalar *array;

 43:   PetscFunctionBegin;
 45:   PetscAssertPointer(file, 2);
 46:   PetscCall(VecGetSize(b, &n));
 47:   PetscCall(VecGetArrayRead(b, &array));
 48:   fprintf(file, "%" PetscInt_FMT "\n", n);
 49:   for (PetscInt i = 0; i < n; i++) fprintf(file, "%" PetscInt_FMT " %16.14e\n", i + 1, array[i]);
 50:   PetscCall(VecRestoreArrayRead(b, &array));
 51:   PetscFunctionReturn(PETSC_SUCCESS);
 52: }