Actual source code: ex249.c

  1: static char help[] = "Test MatCreateSubMatrices\n\n";

  3: #include <petscis.h>
  4: #include <petscmat.h>

  6: int main(int argc, char **args)
  7: {
  8:   Mat         A, *submats, *submats2;
  9:   IS         *irow, *icol;
 10:   PetscInt    i, n;
 11:   PetscMPIInt rank;
 12:   PetscViewer matfd, rowfd, colfd;
 13:   PetscBool   same;
 14:   char        matfile[PETSC_MAX_PATH_LEN], rowfile[PETSC_MAX_PATH_LEN], colfile[PETSC_MAX_PATH_LEN];
 15:   char        rankstr[16] = {0};

 17:   PetscFunctionBeginUser;
 18:   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
 19:   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));

 21:   PetscCall(PetscOptionsGetString(NULL, NULL, "-A", matfile, sizeof(matfile), NULL));
 22:   PetscCall(PetscOptionsGetString(NULL, NULL, "-row", rowfile, sizeof(rowfile), NULL));
 23:   PetscCall(PetscOptionsGetString(NULL, NULL, "-col", colfile, sizeof(colfile), NULL));

 25:   /* Each rank has its own files for row/col ISes */
 26:   PetscCall(PetscSNPrintf(rankstr, 16, "-%d", rank));
 27:   PetscCall(PetscStrlcat(rowfile, rankstr, PETSC_MAX_PATH_LEN));
 28:   PetscCall(PetscStrlcat(colfile, rankstr, PETSC_MAX_PATH_LEN));

 30:   PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD, matfile, FILE_MODE_READ, &matfd));
 31:   PetscCall(PetscViewerBinaryOpen(PETSC_COMM_SELF, rowfile, FILE_MODE_READ, &rowfd));
 32:   PetscCall(PetscViewerBinaryOpen(PETSC_COMM_SELF, colfile, FILE_MODE_READ, &colfd));

 34:   PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
 35:   PetscCall(MatSetFromOptions(A));
 36:   PetscCall(MatLoad(A, matfd));

 38:   /* We stored the number of ISes at the beginning of rowfd */
 39:   PetscCall(PetscViewerBinaryRead(rowfd, &n, 1, NULL, PETSC_INT));
 40:   PetscCall(PetscMalloc2(n, &irow, n, &icol));
 41:   for (i = 0; i < n; i++) {
 42:     PetscCall(ISCreate(PETSC_COMM_SELF, &irow[i]));
 43:     PetscCall(ISCreate(PETSC_COMM_SELF, &icol[i]));
 44:     PetscCall(ISLoad(irow[i], rowfd));
 45:     PetscCall(ISLoad(icol[i], colfd));
 46:   }

 48:   PetscCall(PetscViewerDestroy(&matfd));
 49:   PetscCall(PetscViewerDestroy(&rowfd));
 50:   PetscCall(PetscViewerDestroy(&colfd));

 52:   /* Create submats for the first time */
 53:   PetscCall(MatCreateSubMatrices(A, n, irow, icol, MAT_INITIAL_MATRIX, &submats));

 55:   /* Dup submats to submats2 for later comparison */
 56:   PetscCall(PetscMalloc1(n, &submats2));
 57:   for (i = 0; i < n; i++) PetscCall(MatDuplicate(submats[i], MAT_COPY_VALUES, &submats2[i]));

 59:   /* Create submats again */
 60:   PetscCall(MatCreateSubMatrices(A, n, irow, icol, MAT_REUSE_MATRIX, &submats));

 62:   /* Compare submats and submats2 */
 63:   for (i = 0; i < n; i++) {
 64:     PetscCall(MatEqual(submats[i], submats2[i], &same));
 65:     PetscCheck(same, PETSC_COMM_SELF, PETSC_ERR_PLIB, "submatrix %" PetscInt_FMT " is not same", i);
 66:   }

 68:   PetscCall(MatDestroy(&A));
 69:   for (i = 0; i < n; i++) {
 70:     PetscCall(ISDestroy(&irow[i]));
 71:     PetscCall(ISDestroy(&icol[i]));
 72:   }
 73:   PetscCall(MatDestroySubMatrices(n, &submats));
 74:   PetscCall(MatDestroyMatrices(n, &submats2));
 75:   PetscCall(PetscFree2(irow, icol));
 76:   PetscCall(PetscFinalize());
 77:   return 0;
 78: }

 80: /*TEST

 82:    test:
 83:      suffix: 1
 84:      nsize: 2
 85:      requires: datafilespath double !complex !defined(PETSC_USE_64BIT_INDICES)
 86:      args: -mat_type {{aij baij}} -A ${DATAFILESPATH}/matrices/CreateSubMatrices/A -row ${DATAFILESPATH}/matrices/CreateSubMatrices/row -col ${DATAFILESPATH}/matrices/CreateSubMatrices/col

 88: TEST*/