Actual source code: ex33.c
1: static char help[] = "Tests PetscPartitioner.\n\n";
3: #include <petscpartitioner.h>
5: int main(int argc, char **argv)
6: {
7: PetscPartitioner p;
8: PetscSection partSection, vertexSection = NULL, targetSection = NULL;
9: IS partition, is;
10: PetscMPIInt size, rank;
11: PetscInt nparts, i;
12: PetscInt nv = 4;
13: PetscInt vv[5] = {0, 2, 4, 6, 8};
14: PetscInt vadj[8] = {3, 1, 0, 2, 1, 3, 2, 0};
15: PetscBool sequential;
16: PetscBool vwgts = PETSC_FALSE;
17: PetscBool pwgts = PETSC_FALSE;
19: PetscFunctionBeginUser;
20: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
21: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
22: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
23: nparts = size;
24: PetscCall(PetscOptionsGetInt(NULL, NULL, "-nparts", &nparts, NULL));
25: PetscCall(PetscOptionsGetBool(NULL, NULL, "-vwgts", &vwgts, NULL));
26: PetscCall(PetscOptionsGetBool(NULL, NULL, "-pwgts", &pwgts, NULL));
28: /* create PetscPartitioner */
29: PetscCall(PetscPartitionerCreate(PETSC_COMM_WORLD, &p));
30: PetscCall(PetscPartitionerSetType(p, PETSCPARTITIONERSIMPLE));
31: PetscCall(PetscPartitionerSetFromOptions(p));
33: /* create partition section */
34: PetscCall(PetscSectionCreate(PETSC_COMM_WORLD, &partSection));
36: if (vwgts) { /* create vertex weights section */
37: PetscCall(PetscSectionCreate(PETSC_COMM_WORLD, &vertexSection));
38: PetscCall(PetscSectionSetChart(vertexSection, 0, nv));
39: for (i = 0; i < nv; i++) PetscCall(PetscSectionSetDof(vertexSection, i, 1));
40: PetscCall(PetscSectionSetUp(vertexSection));
41: }
43: if (pwgts) { /* create partition weights section */
44: PetscCall(PetscSectionCreate(PETSC_COMM_WORLD, &targetSection));
45: PetscCall(PetscSectionSetChart(targetSection, 0, nparts));
46: for (i = 0; i < nparts; i++) PetscCall(PetscSectionSetDof(targetSection, i, 1));
47: PetscCall(PetscSectionSetUp(targetSection));
48: }
50: if (PetscDefined(USE_LOG)) { /* Test logging */
51: PetscLogEvent event;
53: PetscCall(PetscLogEventRegister("MyPartitionerEvent", PETSCPARTITIONER_CLASSID, &event));
54: { /* PetscLogEventExcludeClass is broken, new events are not deactivated */
55: char logList[256];
56: PetscBool opt, pkg;
58: PetscCall(PetscOptionsGetString(NULL, NULL, "-log_exclude", logList, sizeof(logList), &opt));
59: if (opt) {
60: PetscCall(PetscStrInList("partitioner", logList, ',', &pkg));
61: if (pkg) PetscCall(PetscLogEventExcludeClass(PETSCPARTITIONER_CLASSID));
62: }
63: }
64: PetscCall(PetscLogEventBegin(event, p, NULL, NULL, NULL));
65: PetscCall(PetscLogEventEnd(event, p, NULL, NULL, NULL));
66: }
68: /* test setup and reset */
69: PetscCall(PetscPartitionerSetUp(p));
70: PetscCall(PetscPartitionerReset(p));
72: /* test partitioning an empty graph */
73: PetscCall(PetscPartitionerPartition(p, nparts, 0, NULL, NULL, vertexSection, NULL, targetSection, partSection, &partition));
74: PetscCall(PetscObjectSetName((PetscObject)partSection, "NULL SECTION"));
75: PetscCall(PetscSectionView(partSection, NULL));
76: PetscCall(ISOnComm(partition, PETSC_COMM_WORLD, PETSC_USE_POINTER, &is));
77: PetscCall(PetscObjectSetName((PetscObject)is, "NULL PARTITION"));
78: PetscCall(ISView(is, NULL));
79: PetscCall(ISDestroy(&is));
80: PetscCall(ISDestroy(&partition));
82: /* test view from options */
83: PetscCall(PetscPartitionerViewFromOptions(p, NULL, "-part_view"));
85: /* test partitioning a graph on one process only (not main) */
86: if (rank == size - 1) {
87: PetscCall(PetscPartitionerPartition(p, nparts, nv, vv, vadj, vertexSection, NULL, targetSection, partSection, &partition));
88: } else {
89: PetscCall(PetscPartitionerPartition(p, nparts, 0, NULL, NULL, vertexSection, NULL, targetSection, partSection, &partition));
90: }
91: PetscCall(PetscObjectSetName((PetscObject)partSection, "SEQ SECTION"));
92: PetscCall(PetscSectionView(partSection, NULL));
93: PetscCall(ISOnComm(partition, PETSC_COMM_WORLD, PETSC_USE_POINTER, &is));
94: PetscCall(PetscObjectSetName((PetscObject)is, "SEQ PARTITION"));
95: PetscCall(ISView(is, NULL));
96: PetscCall(ISDestroy(&is));
97: PetscCall(ISDestroy(&partition));
99: PetscCall(PetscObjectTypeCompareAny((PetscObject)p, &sequential, PETSCPARTITIONERCHACO, NULL));
100: if (sequential) goto finally;
102: /* test partitioning a graph on a subset of the processes only */
103: if (rank % 2) {
104: PetscCall(PetscPartitionerPartition(p, nparts, 0, NULL, NULL, NULL, NULL, targetSection, partSection, &partition));
105: } else {
106: PetscInt i, totv = nv * ((size + 1) / 2), *pvadj;
108: PetscCall(PetscMalloc1(2 * nv, &pvadj));
109: for (i = 0; i < nv; i++) {
110: pvadj[2 * i] = (nv * (rank / 2) + totv + i - 1) % totv;
111: pvadj[2 * i + 1] = (nv * (rank / 2) + totv + i + 1) % totv;
112: }
113: PetscCall(PetscPartitionerPartition(p, nparts, nv, vv, pvadj, NULL, NULL, targetSection, partSection, &partition));
114: PetscCall(PetscFree(pvadj));
115: }
116: PetscCall(PetscObjectSetName((PetscObject)partSection, "PARVOID SECTION"));
117: PetscCall(PetscSectionView(partSection, NULL));
118: PetscCall(ISOnComm(partition, PETSC_COMM_WORLD, PETSC_USE_POINTER, &is));
119: PetscCall(PetscObjectSetName((PetscObject)is, "PARVOID PARTITION"));
120: PetscCall(ISView(is, NULL));
121: PetscCall(ISDestroy(&is));
122: PetscCall(ISDestroy(&partition));
124: finally:
125: PetscCall(PetscSectionDestroy(&partSection));
126: PetscCall(PetscSectionDestroy(&vertexSection));
127: PetscCall(PetscSectionDestroy(&targetSection));
128: PetscCall(PetscPartitionerDestroy(&p));
129: PetscCall(PetscFinalize());
130: return 0;
131: }
133: /*TEST
135: test:
136: suffix: default
138: testset:
139: requires: defined(PETSC_USE_LOG)
140: args: -petscpartitioner_type simple -log_view
141: filter: grep MyPartitionerEvent | cut -d " " -f 1
142: test:
143: suffix: log_include
144: test:
145: suffix: log_exclude
146: args: -log_exclude partitioner
147: output_file: output/empty.out
149: test:
150: suffix: simple
151: nsize: {{1 2 3}separate output}
152: args: -nparts {{1 2 3}separate output} -pwgts {{false true}separate output} -petscpartitioner_type simple -petscpartitioner_view
154: test:
155: suffix: shell
156: nsize: {{1 2 3}separate output}
157: args: -nparts {{1 2 3}separate output} -petscpartitioner_type shell -petscpartitioner_shell_random -petscpartitioner_view
159: test:
160: suffix: gather
161: nsize: {{1 2 3}separate output}
162: args: -nparts {{1 2 3}separate output} -petscpartitioner_type gather -petscpartitioner_view -petscpartitioner_view_graph
164: test:
165: requires: parmetis
166: suffix: parmetis
167: nsize: {{1 2 3}separate output}
168: args: -nparts {{1 2 3}separate output} -pwgts {{false true}} -vwgts {{false true}}
169: args: -petscpartitioner_type parmetis -petscpartitioner_view -petscpartitioner_view_graph
171: test:
172: requires: parmetis
173: suffix: parmetis_type
174: nsize: {{1 2}}
175: args: -petscpartitioner_type parmetis -part_view
176: args: -petscpartitioner_parmetis_type {{kway rb}separate output}
177: filter: grep "ParMetis type"
179: test:
180: requires: ptscotch
181: suffix: ptscotch
182: nsize: {{1 2 3}separate output}
183: args: -nparts {{1 2 3}separate output} -pwgts {{false true}separate output} -vwgts {{false true}}
184: args: -petscpartitioner_type ptscotch -petscpartitioner_view -petscpartitioner_view_graph
186: test:
187: requires: ptscotch
188: suffix: ptscotch_strategy
189: nsize: {{1 2}}
190: args: -petscpartitioner_type ptscotch -part_view
191: args: -petscpartitioner_ptscotch_strategy {{DEFAULT QUALITY SPEED BALANCE SAFETY SCALABILITY RECURSIVE REMAP}separate output}
192: filter: grep "partitioning strategy"
194: test:
195: requires: chaco
196: suffix: chaco
197: nsize: {{1 2 3}separate output}
198: args: -nparts {{1}separate output} -petscpartitioner_type chaco -petscpartitioner_view -petscpartitioner_view_graph
200: test:
201: TODO: non reproducible (uses C stdlib rand())
202: requires: chaco
203: suffix: chaco
204: nsize: {{1 2 3}separate output}
205: args: -nparts {{2 3}separate output} -petscpartitioner_type chaco -petscpartitioner_view -petscpartitioner_view_graph
207: TEST*/