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, NULL, 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*/