Actual source code: ex87.c


  2: static char help[] = "Tests MatCreateSubMatrices() for SBAIJ matrices\n\n";

  4: #include <petscmat.h>

  6: int main(int argc, char **args)
  7: {
  8:   Mat         BAIJ, SBAIJ, *subBAIJ, *subSBAIJ;
  9:   PetscViewer viewer;
 10:   char        file[PETSC_MAX_PATH_LEN];
 11:   PetscBool   flg;
 12:   PetscInt    n = 2, issize, M, N;
 13:   PetscMPIInt rank;
 14:   IS          isrow, iscol, irow[n], icol[n];

 16:   PetscFunctionBeginUser;
 17:   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
 18:   PetscCall(PetscOptionsGetString(NULL, NULL, "-f", file, sizeof(file), &flg));
 19:   PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &viewer));

 21:   PetscCall(MatCreate(PETSC_COMM_WORLD, &BAIJ));
 22:   PetscCall(MatSetType(BAIJ, MATMPIBAIJ));
 23:   PetscCall(MatLoad(BAIJ, viewer));
 24:   PetscCall(PetscViewerDestroy(&viewer));

 26:   PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, file, FILE_MODE_READ, &viewer));
 27:   PetscCall(MatCreate(PETSC_COMM_WORLD, &SBAIJ));
 28:   PetscCall(MatSetType(SBAIJ, MATMPISBAIJ));
 29:   PetscCall(MatLoad(SBAIJ, viewer));
 30:   PetscCall(PetscViewerDestroy(&viewer));

 32:   PetscCall(MatGetSize(BAIJ, &M, &N));
 33:   issize = M / 4;
 34:   PetscCall(ISCreateStride(PETSC_COMM_SELF, issize, 0, 1, &isrow));
 35:   irow[0] = irow[1] = isrow;
 36:   issize            = N / 2;
 37:   PetscCall(ISCreateStride(PETSC_COMM_SELF, issize, 0, 1, &iscol));
 38:   icol[0] = icol[1] = iscol;
 39:   PetscCall(MatCreateSubMatrices(BAIJ, n, irow, icol, MAT_INITIAL_MATRIX, &subBAIJ));
 40:   PetscCall(MatCreateSubMatrices(BAIJ, n, irow, icol, MAT_REUSE_MATRIX, &subBAIJ));

 42:   /* irow and icol must be same for SBAIJ matrices! */
 43:   icol[0] = icol[1] = isrow;
 44:   PetscCall(MatCreateSubMatrices(SBAIJ, n, irow, icol, MAT_INITIAL_MATRIX, &subSBAIJ));
 45:   PetscCall(MatCreateSubMatrices(SBAIJ, n, irow, icol, MAT_REUSE_MATRIX, &subSBAIJ));

 47:   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
 48:   if (rank == 0) {
 49:     PetscCall(MatView(subBAIJ[0], PETSC_VIEWER_STDOUT_SELF));
 50:     PetscCall(MatView(subSBAIJ[0], PETSC_VIEWER_STDOUT_SELF));
 51:   }

 53:   /* Free data structures */
 54:   PetscCall(ISDestroy(&isrow));
 55:   PetscCall(ISDestroy(&iscol));
 56:   PetscCall(MatDestroySubMatrices(n, &subBAIJ));
 57:   PetscCall(MatDestroySubMatrices(n, &subSBAIJ));
 58:   PetscCall(MatDestroy(&BAIJ));
 59:   PetscCall(MatDestroy(&SBAIJ));

 61:   PetscCall(PetscFinalize());
 62:   return 0;
 63: }