1 static char help[] = "Tests PetscPartitioner.\n\n"; 2 3 #include <petscpartitioner.h> 4 5 int main(int argc, char **argv) { 6 PetscPartitioner p; 7 PetscSection partSection, vertexSection = NULL, targetSection = NULL; 8 IS partition, is; 9 PetscMPIInt size, rank; 10 PetscInt nparts, i; 11 PetscInt nv = 4; 12 PetscInt vv[5] = {0, 2, 4, 6, 8}; 13 PetscInt vadj[8] = {3, 1, 0, 2, 1, 3, 2, 0}; 14 PetscBool sequential; 15 PetscBool vwgts = PETSC_FALSE; 16 PetscBool pwgts = PETSC_FALSE; 17 18 PetscFunctionBeginUser; 19 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 20 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 21 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 22 nparts = size; 23 PetscCall(PetscOptionsGetInt(NULL, NULL, "-nparts", &nparts, NULL)); 24 PetscCall(PetscOptionsGetBool(NULL, NULL, "-vwgts", &vwgts, NULL)); 25 PetscCall(PetscOptionsGetBool(NULL, NULL, "-pwgts", &pwgts, NULL)); 26 27 /* create PetscPartitioner */ 28 PetscCall(PetscPartitionerCreate(PETSC_COMM_WORLD, &p)); 29 PetscCall(PetscPartitionerSetType(p, PETSCPARTITIONERSIMPLE)); 30 PetscCall(PetscPartitionerSetFromOptions(p)); 31 32 /* create partition section */ 33 PetscCall(PetscSectionCreate(PETSC_COMM_WORLD, &partSection)); 34 35 if (vwgts) { /* create vertex weights section */ 36 PetscCall(PetscSectionCreate(PETSC_COMM_WORLD, &vertexSection)); 37 PetscCall(PetscSectionSetChart(vertexSection, 0, nv)); 38 for (i = 0; i < nv; i++) PetscCall(PetscSectionSetDof(vertexSection, i, 1)); 39 PetscCall(PetscSectionSetUp(vertexSection)); 40 } 41 42 if (pwgts) { /* create partition weights section */ 43 PetscCall(PetscSectionCreate(PETSC_COMM_WORLD, &targetSection)); 44 PetscCall(PetscSectionSetChart(targetSection, 0, nparts)); 45 for (i = 0; i < nparts; i++) PetscCall(PetscSectionSetDof(targetSection, i, 1)); 46 PetscCall(PetscSectionSetUp(targetSection)); 47 } 48 49 #if defined(PETSC_USE_LOG) 50 { /* Test logging */ 51 PetscLogEvent event; 52 53 PetscCall(PetscLogEventRegister("MyPartitionerEvent", PETSCPARTITIONER_CLASSID, &event)); 54 { /* PetscLogEventExcludeClass is broken, new events are not deactivated */ 55 char logList[256]; 56 PetscBool opt, pkg; 57 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 } 67 #endif 68 69 /* test setup and reset */ 70 PetscCall(PetscPartitionerSetUp(p)); 71 PetscCall(PetscPartitionerReset(p)); 72 73 /* test partitioning an empty graph */ 74 PetscCall(PetscPartitionerPartition(p, nparts, 0, NULL, NULL, vertexSection, targetSection, partSection, &partition)); 75 PetscCall(PetscObjectSetName((PetscObject)partSection, "NULL SECTION")); 76 PetscCall(PetscSectionView(partSection, NULL)); 77 PetscCall(ISOnComm(partition, PETSC_COMM_WORLD, PETSC_USE_POINTER, &is)); 78 PetscCall(PetscObjectSetName((PetscObject)is, "NULL PARTITION")); 79 PetscCall(ISView(is, NULL)); 80 PetscCall(ISDestroy(&is)); 81 PetscCall(ISDestroy(&partition)); 82 83 /* test view from options */ 84 PetscCall(PetscPartitionerViewFromOptions(p, NULL, "-part_view")); 85 86 /* test partitioning a graph on one process only (not main) */ 87 if (rank == size - 1) { 88 PetscCall(PetscPartitionerPartition(p, nparts, nv, vv, vadj, vertexSection, targetSection, partSection, &partition)); 89 } else { 90 PetscCall(PetscPartitionerPartition(p, nparts, 0, NULL, NULL, vertexSection, targetSection, partSection, &partition)); 91 } 92 PetscCall(PetscObjectSetName((PetscObject)partSection, "SEQ SECTION")); 93 PetscCall(PetscSectionView(partSection, NULL)); 94 PetscCall(ISOnComm(partition, PETSC_COMM_WORLD, PETSC_USE_POINTER, &is)); 95 PetscCall(PetscObjectSetName((PetscObject)is, "SEQ PARTITION")); 96 PetscCall(ISView(is, NULL)); 97 PetscCall(ISDestroy(&is)); 98 PetscCall(ISDestroy(&partition)); 99 100 PetscCall(PetscObjectTypeCompareAny((PetscObject)p, &sequential, PETSCPARTITIONERCHACO, NULL)); 101 if (sequential) goto finally; 102 103 /* test partitioning a graph on a subset of the processess only */ 104 if (rank % 2) { 105 PetscCall(PetscPartitionerPartition(p, nparts, 0, NULL, NULL, NULL, targetSection, partSection, &partition)); 106 } else { 107 PetscInt i, totv = nv * ((size + 1) / 2), *pvadj; 108 109 PetscCall(PetscMalloc1(2 * nv, &pvadj)); 110 for (i = 0; i < nv; i++) { 111 pvadj[2 * i] = (nv * (rank / 2) + totv + i - 1) % totv; 112 pvadj[2 * i + 1] = (nv * (rank / 2) + totv + i + 1) % totv; 113 } 114 PetscCall(PetscPartitionerPartition(p, nparts, nv, vv, pvadj, NULL, targetSection, partSection, &partition)); 115 PetscCall(PetscFree(pvadj)); 116 } 117 PetscCall(PetscObjectSetName((PetscObject)partSection, "PARVOID SECTION")); 118 PetscCall(PetscSectionView(partSection, NULL)); 119 PetscCall(ISOnComm(partition, PETSC_COMM_WORLD, PETSC_USE_POINTER, &is)); 120 PetscCall(PetscObjectSetName((PetscObject)is, "PARVOID PARTITION")); 121 PetscCall(ISView(is, NULL)); 122 PetscCall(ISDestroy(&is)); 123 PetscCall(ISDestroy(&partition)); 124 125 finally: 126 PetscCall(PetscSectionDestroy(&partSection)); 127 PetscCall(PetscSectionDestroy(&vertexSection)); 128 PetscCall(PetscSectionDestroy(&targetSection)); 129 PetscCall(PetscPartitionerDestroy(&p)); 130 PetscCall(PetscFinalize()); 131 return 0; 132 } 133 134 /*TEST 135 136 test: 137 suffix: default 138 139 testset: 140 requires: defined(PETSC_USE_LOG) 141 args: -petscpartitioner_type simple -log_summary 142 filter: grep MyPartitionerEvent | cut -d " " -f 1 143 test: 144 suffix: log_include 145 test: 146 suffix: log_exclude 147 args: -log_exclude partitioner 148 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 153 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 158 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 163 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 170 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" 178 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 185 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" 193 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 199 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 206 207 TEST*/ 208