1 static char help[] = "Test MatCreateSubMatrices\n\n"; 2 3 #include <petscis.h> 4 #include <petscmat.h> 5 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}; 16 17 PetscFunctionBeginUser; 18 PetscCall(PetscInitialize(&argc, &args, NULL, help)); 19 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 20 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)); 24 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)); 29 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)); 33 34 PetscCall(MatCreate(PETSC_COMM_WORLD, &A)); 35 PetscCall(MatSetFromOptions(A)); 36 PetscCall(MatLoad(A, matfd)); 37 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 } 47 48 PetscCall(PetscViewerDestroy(&matfd)); 49 PetscCall(PetscViewerDestroy(&rowfd)); 50 PetscCall(PetscViewerDestroy(&colfd)); 51 52 /* Create submats for the first time */ 53 PetscCall(MatCreateSubMatrices(A, n, irow, icol, MAT_INITIAL_MATRIX, &submats)); 54 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])); 58 59 /* Create submats again */ 60 PetscCall(MatCreateSubMatrices(A, n, irow, icol, MAT_REUSE_MATRIX, &submats)); 61 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 } 67 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 } 79 80 /*TEST 81 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 87 output_file: output/empty.out 88 89 TEST*/ 90