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 */
12
13 #include <petscmat.h>
14
main(int argc,char ** argv)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;
36
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));
46
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));
50
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));
56
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));
63
64 PetscCall(MatPermute(A, isrow, iscol, &B));
65 PetscCall(PetscViewerASCIIPrintf(viewer, "Permuted matrix\n"));
66 PetscCall(MatView(B, viewer));
67
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));
73
74 /* Free data structures */
75 PetscCall(ISDestroy(&isrow));
76 PetscCall(ISDestroy(&iscol));
77 PetscCall(MatDestroy(&A));
78 PetscCall(MatDestroy(&B));
79
80 PetscCall(PetscFinalize());
81 return 0;
82 }
83
84 /*TEST
85
86 build:
87 requires: !complex
88
89 test:
90 args: -view_sparse
91
92 test:
93 suffix: 2
94 nsize: 2
95 args: -view_sparse
96
97 test:
98 suffix: 2b
99 nsize: 2
100 args: -mat_type baij -view_sparse
101
102 test:
103 suffix: 3
104 nsize: 3
105 args: -view_sparse
106
107 test:
108 suffix: 3b
109 nsize: 3
110 args: -mat_type baij -view_sparse
111
112 test:
113 suffix: dense
114 args: -mat_type dense
115
116 TEST*/
117