Actual source code: ex166.c
1: static char help[] = "Tests MatPermute() for a square matrix in parallel.\n\n";
2: #include <petscmat.h>
3: /* Results:
4: Sequential:
5: - seqaij: correct permutation
6: - seqbaij: permutation not supported for this MATTYPE
7: - seqsbaij: permutation not supported for this MATTYPE
8: Parallel:
9: - mpiaij: incorrect permutation (disable this op for now)
10: - seqbaij: correct permutation, but all manner of memory leaks
11: (disable this op for now, because it appears broken for nonsquare matrices; see ex151)
12: - seqsbaij: permutation not supported for this MATTYPE
14: */
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, 4, 4.},
25: {3, 0, 5.},
26: {3, 3, 6.},
27: {4, 1, 7.},
28: {4, 4, 8.}
29: };
30: const PetscInt ixrow[5] = {4, 2, 1, 3, 0}, ixcol[5] = {3, 2, 1, 4, 0};
31: Mat A, B;
32: PetscInt i, rstart, rend, cstart, cend;
33: IS isrow, iscol;
34: PetscViewer viewer, sviewer;
35: PetscBool view_sparse;
38: PetscInitialize(&argc, &argv, (char *)0, help);
39: /* ------- Assemble matrix, --------- */
40: MatCreate(PETSC_COMM_WORLD, &A);
41: MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, 5, 5);
42: MatSetFromOptions(A);
43: MatSetUp(A);
44: MatGetOwnershipRange(A, &rstart, &rend);
45: MatGetOwnershipRangeColumn(A, &cstart, &cend);
47: for (i = 0; i < (PetscInt)PETSC_STATIC_ARRAY_LENGTH(entries); i++) MatSetValue(A, entries[i].i, entries[i].j, entries[i].v, INSERT_VALUES);
48: MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
49: MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
51: /* ------ Prepare index sets ------ */
52: ISCreateGeneral(PETSC_COMM_WORLD, rend - rstart, ixrow + rstart, PETSC_USE_POINTER, &isrow);
53: ISCreateGeneral(PETSC_COMM_SELF, 5, ixcol, PETSC_USE_POINTER, &iscol);
54: ISSetPermutation(isrow);
55: ISSetPermutation(iscol);
57: PetscViewerASCIIGetStdout(PETSC_COMM_WORLD, &viewer);
58: view_sparse = PETSC_FALSE;
59: PetscOptionsGetBool(NULL, NULL, "-view_sparse", &view_sparse, NULL);
60: if (!view_sparse) PetscViewerPushFormat(viewer, PETSC_VIEWER_ASCII_DENSE);
61: PetscViewerASCIIPrintf(viewer, "Original matrix\n");
62: MatView(A, viewer);
64: MatPermute(A, isrow, iscol, &B);
65: PetscViewerASCIIPrintf(viewer, "Permuted matrix\n");
66: MatView(B, viewer);
67: if (!view_sparse) PetscViewerPopFormat(viewer);
69: PetscViewerASCIIPrintf(viewer, "Row permutation\n");
70: ISView(isrow, viewer);
71: PetscViewerASCIIPrintf(viewer, "Column permutation\n");
72: PetscViewerGetSubViewer(viewer, PETSC_COMM_SELF, &sviewer);
73: ISView(iscol, sviewer);
74: PetscViewerRestoreSubViewer(viewer, PETSC_COMM_SELF, &sviewer);
76: /* Free data structures */
77: ISDestroy(&isrow);
78: ISDestroy(&iscol);
79: MatDestroy(&A);
80: MatDestroy(&B);
82: PetscFinalize();
83: return 0;
84: }