1 #include <petscmat.h>
2 #include <petsc/private/petscimpl.h>
3
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.
9
10 */
MatDumpSPAI(Mat A,FILE * file)11 PetscErrorCode MatDumpSPAI(Mat A, FILE *file)
12 {
13 PetscMPIInt size;
14 PetscInt n;
15 MPI_Comm comm;
16
17 PetscFunctionBegin;
18 PetscValidHeaderSpecific(A, MAT_CLASSID, 1);
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;
30
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 }
37
VecDumpSPAI(Vec b,FILE * file)38 PetscErrorCode VecDumpSPAI(Vec b, FILE *file)
39 {
40 PetscInt n;
41 const PetscScalar *array;
42
43 PetscFunctionBegin;
44 PetscValidHeaderSpecific(b, VEC_CLASSID, 1);
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 }
53