Actual source code: ex151.c

  1: static char help[] = "Tests MatPermute() in parallel.\n\n";
  2: /* Results:
  3:    Sequential:
  4:    - seqaij:   correct permutation
  5:    - seqbaij:  permutation not supported for this MATTYPE
  6:    - seqsbaij: permutation not supported for this MATTYPE
  7:    Parallel:
  8:    - mpiaij:   correct permutation
  9:    - mpibaij:  correct permutation
 10:    - mpisbaij: permutation not supported for this MATTYPE
 11:  */

 13: #include <petscmat.h>

 15: int main(int argc, char **argv)
 16: {
 17:   const struct {
 18:     PetscInt    i, j;
 19:     PetscScalar v;
 20:   } entries[] = {
 21:     {0, 3, 1.},
 22:     {1, 2, 2.},
 23:     {2, 1, 3.},
 24:     {2, 5, 4.},
 25:     {3, 0, 5.},
 26:     {3, 6, 6.},
 27:     {4, 1, 7.},
 28:     {4, 4, 8.}
 29:   };
 30:   const PetscInt ixrow[5] = {4, 2, 1, 0, 3}, ixcol[7] = {5, 3, 6, 1, 2, 0, 4};
 31:   Mat            A, B;
 32:   PetscInt       i, rstart, rend, cstart, cend;
 33:   IS             isrow, iscol;
 34:   PetscViewer    viewer;
 35:   PetscBool      view_sparse;

 37:   PetscFunctionBeginUser;
 38:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 39:   /* ------- Assemble matrix, --------- */
 40:   PetscCall(MatCreate(PETSC_COMM_WORLD, &A));
 41:   PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, 5, 7));
 42:   PetscCall(MatSetFromOptions(A));
 43:   PetscCall(MatSetUp(A));
 44:   PetscCall(MatGetOwnershipRange(A, &rstart, &rend));
 45:   PetscCall(MatGetOwnershipRangeColumn(A, &cstart, &cend));

 47:   for (i = 0; i < (PetscInt)PETSC_STATIC_ARRAY_LENGTH(entries); i++) PetscCall(MatSetValue(A, entries[i].i, entries[i].j, entries[i].v, INSERT_VALUES));
 48:   PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
 49:   PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));

 51:   /* ------ Prepare index sets ------ */
 52:   PetscCall(ISCreateGeneral(PETSC_COMM_WORLD, rend - rstart, ixrow + rstart, PETSC_USE_POINTER, &isrow));
 53:   PetscCall(ISCreateGeneral(PETSC_COMM_WORLD, cend - cstart, ixcol + cstart, PETSC_USE_POINTER, &iscol));
 54:   PetscCall(ISSetPermutation(isrow));
 55:   PetscCall(ISSetPermutation(iscol));

 57:   PetscCall(PetscViewerASCIIGetStdout(PETSC_COMM_WORLD, &viewer));
 58:   view_sparse = PETSC_FALSE;
 59:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-view_sparse", &view_sparse, NULL));
 60:   if (!view_sparse) PetscCall(PetscViewerPushFormat(viewer, PETSC_VIEWER_ASCII_DENSE));
 61:   PetscCall(PetscViewerASCIIPrintf(viewer, "Original matrix\n"));
 62:   PetscCall(MatView(A, viewer));

 64:   PetscCall(MatPermute(A, isrow, iscol, &B));
 65:   PetscCall(PetscViewerASCIIPrintf(viewer, "Permuted matrix\n"));
 66:   PetscCall(MatView(B, viewer));

 68:   if (!view_sparse) PetscCall(PetscViewerPopFormat(viewer));
 69:   PetscCall(PetscViewerASCIIPrintf(viewer, "Row permutation\n"));
 70:   PetscCall(ISView(isrow, viewer));
 71:   PetscCall(PetscViewerASCIIPrintf(viewer, "Column permutation\n"));
 72:   PetscCall(ISView(iscol, viewer));

 74:   /* Free data structures */
 75:   PetscCall(ISDestroy(&isrow));
 76:   PetscCall(ISDestroy(&iscol));
 77:   PetscCall(MatDestroy(&A));
 78:   PetscCall(MatDestroy(&B));

 80:   PetscCall(PetscFinalize());
 81:   return 0;
 82: }

 84: /*TEST

 86:    build:
 87:       requires: !complex

 89:    test:
 90:       args: -view_sparse

 92:    test:
 93:       suffix: 2
 94:       nsize: 2
 95:       args: -view_sparse

 97:    test:
 98:       suffix: 2b
 99:       nsize: 2
100:       args: -mat_type baij -view_sparse

102:    test:
103:       suffix: 3
104:       nsize: 3
105:       args: -view_sparse

107:    test:
108:       suffix: 3b
109:       nsize: 3
110:       args: -mat_type baij -view_sparse

112:    test:
113:       suffix: dense
114:       args: -mat_type dense

116: TEST*/