1 static char help[] = "Tests PetscPartitioner.\n\n"; 2 3 #include <petscpartitioner.h> 4 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; 18 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)); 27 28 /* create PetscPartitioner */ 29 PetscCall(PetscPartitionerCreate(PETSC_COMM_WORLD, &p)); 30 PetscCall(PetscPartitionerSetType(p, PETSCPARTITIONERSIMPLE)); 31 PetscCall(PetscPartitionerSetFromOptions(p)); 32 33 /* create partition section */ 34 PetscCall(PetscSectionCreate(PETSC_COMM_WORLD, &partSection)); 35 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 } 42 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 } 49 50 #if defined(PETSC_USE_LOG) 51 { /* Test logging */ 52 PetscLogEvent event; 53 54 PetscCall(PetscLogEventRegister("MyPartitionerEvent", PETSCPARTITIONER_CLASSID, &event)); 55 { /* PetscLogEventExcludeClass is broken, new events are not deactivated */ 56 char logList[256]; 57 PetscBool opt, pkg; 58 59 PetscCall(PetscOptionsGetString(NULL, NULL, "-log_exclude", logList, sizeof(logList), &opt)); 60 if (opt) { 61 PetscCall(PetscStrInList("partitioner", logList, ',', &pkg)); 62 if (pkg) PetscCall(PetscLogEventExcludeClass(PETSCPARTITIONER_CLASSID)); 63 } 64 } 65 PetscCall(PetscLogEventBegin(event, p, NULL, NULL, NULL)); 66 PetscCall(PetscLogEventEnd(event, p, NULL, NULL, NULL)); 67 } 68 #endif 69 70 /* test setup and reset */ 71 PetscCall(PetscPartitionerSetUp(p)); 72 PetscCall(PetscPartitionerReset(p)); 73 74 /* test partitioning an empty graph */ 75 PetscCall(PetscPartitionerPartition(p, nparts, 0, NULL, NULL, vertexSection, targetSection, partSection, &partition)); 76 PetscCall(PetscObjectSetName((PetscObject)partSection, "NULL SECTION")); 77 PetscCall(PetscSectionView(partSection, NULL)); 78 PetscCall(ISOnComm(partition, PETSC_COMM_WORLD, PETSC_USE_POINTER, &is)); 79 PetscCall(PetscObjectSetName((PetscObject)is, "NULL PARTITION")); 80 PetscCall(ISView(is, NULL)); 81 PetscCall(ISDestroy(&is)); 82 PetscCall(ISDestroy(&partition)); 83 84 /* test view from options */ 85 PetscCall(PetscPartitionerViewFromOptions(p, NULL, "-part_view")); 86 87 /* test partitioning a graph on one process only (not main) */ 88 if (rank == size - 1) { 89 PetscCall(PetscPartitionerPartition(p, nparts, nv, vv, vadj, vertexSection, targetSection, partSection, &partition)); 90 } else { 91 PetscCall(PetscPartitionerPartition(p, nparts, 0, NULL, NULL, vertexSection, targetSection, partSection, &partition)); 92 } 93 PetscCall(PetscObjectSetName((PetscObject)partSection, "SEQ SECTION")); 94 PetscCall(PetscSectionView(partSection, NULL)); 95 PetscCall(ISOnComm(partition, PETSC_COMM_WORLD, PETSC_USE_POINTER, &is)); 96 PetscCall(PetscObjectSetName((PetscObject)is, "SEQ PARTITION")); 97 PetscCall(ISView(is, NULL)); 98 PetscCall(ISDestroy(&is)); 99 PetscCall(ISDestroy(&partition)); 100 101 PetscCall(PetscObjectTypeCompareAny((PetscObject)p, &sequential, PETSCPARTITIONERCHACO, NULL)); 102 if (sequential) goto finally; 103 104 /* test partitioning a graph on a subset of the processes only */ 105 if (rank % 2) { 106 PetscCall(PetscPartitionerPartition(p, nparts, 0, NULL, NULL, NULL, targetSection, partSection, &partition)); 107 } else { 108 PetscInt i, totv = nv * ((size + 1) / 2), *pvadj; 109 110 PetscCall(PetscMalloc1(2 * nv, &pvadj)); 111 for (i = 0; i < nv; i++) { 112 pvadj[2 * i] = (nv * (rank / 2) + totv + i - 1) % totv; 113 pvadj[2 * i + 1] = (nv * (rank / 2) + totv + i + 1) % totv; 114 } 115 PetscCall(PetscPartitionerPartition(p, nparts, nv, vv, pvadj, NULL, targetSection, partSection, &partition)); 116 PetscCall(PetscFree(pvadj)); 117 } 118 PetscCall(PetscObjectSetName((PetscObject)partSection, "PARVOID SECTION")); 119 PetscCall(PetscSectionView(partSection, NULL)); 120 PetscCall(ISOnComm(partition, PETSC_COMM_WORLD, PETSC_USE_POINTER, &is)); 121 PetscCall(PetscObjectSetName((PetscObject)is, "PARVOID PARTITION")); 122 PetscCall(ISView(is, NULL)); 123 PetscCall(ISDestroy(&is)); 124 PetscCall(ISDestroy(&partition)); 125 126 finally: 127 PetscCall(PetscSectionDestroy(&partSection)); 128 PetscCall(PetscSectionDestroy(&vertexSection)); 129 PetscCall(PetscSectionDestroy(&targetSection)); 130 PetscCall(PetscPartitionerDestroy(&p)); 131 PetscCall(PetscFinalize()); 132 return 0; 133 } 134 135 /*TEST 136 137 test: 138 suffix: default 139 140 testset: 141 requires: defined(PETSC_USE_LOG) 142 args: -petscpartitioner_type simple -log_summary 143 filter: grep MyPartitionerEvent | cut -d " " -f 1 144 test: 145 suffix: log_include 146 test: 147 suffix: log_exclude 148 args: -log_exclude partitioner 149 150 test: 151 suffix: simple 152 nsize: {{1 2 3}separate output} 153 args: -nparts {{1 2 3}separate output} -pwgts {{false true}separate output} -petscpartitioner_type simple -petscpartitioner_view 154 155 test: 156 suffix: shell 157 nsize: {{1 2 3}separate output} 158 args: -nparts {{1 2 3}separate output} -petscpartitioner_type shell -petscpartitioner_shell_random -petscpartitioner_view 159 160 test: 161 suffix: gather 162 nsize: {{1 2 3}separate output} 163 args: -nparts {{1 2 3}separate output} -petscpartitioner_type gather -petscpartitioner_view -petscpartitioner_view_graph 164 165 test: 166 requires: parmetis 167 suffix: parmetis 168 nsize: {{1 2 3}separate output} 169 args: -nparts {{1 2 3}separate output} -pwgts {{false true}} -vwgts {{false true}} 170 args: -petscpartitioner_type parmetis -petscpartitioner_view -petscpartitioner_view_graph 171 172 test: 173 requires: parmetis 174 suffix: parmetis_type 175 nsize: {{1 2}} 176 args: -petscpartitioner_type parmetis -part_view 177 args: -petscpartitioner_parmetis_type {{kway rb}separate output} 178 filter: grep "ParMetis type" 179 180 test: 181 requires: ptscotch 182 suffix: ptscotch 183 nsize: {{1 2 3}separate output} 184 args: -nparts {{1 2 3}separate output} -pwgts {{false true}separate output} -vwgts {{false true}} 185 args: -petscpartitioner_type ptscotch -petscpartitioner_view -petscpartitioner_view_graph 186 187 test: 188 requires: ptscotch 189 suffix: ptscotch_strategy 190 nsize: {{1 2}} 191 args: -petscpartitioner_type ptscotch -part_view 192 args: -petscpartitioner_ptscotch_strategy {{DEFAULT QUALITY SPEED BALANCE SAFETY SCALABILITY RECURSIVE REMAP}separate output} 193 filter: grep "partitioning strategy" 194 195 test: 196 requires: chaco 197 suffix: chaco 198 nsize: {{1 2 3}separate output} 199 args: -nparts {{1}separate output} -petscpartitioner_type chaco -petscpartitioner_view -petscpartitioner_view_graph 200 201 test: 202 TODO: non reproducible (uses C stdlib rand()) 203 requires: chaco 204 suffix: chaco 205 nsize: {{1 2 3}separate output} 206 args: -nparts {{2 3}separate output} -petscpartitioner_type chaco -petscpartitioner_view -petscpartitioner_view_graph 207 208 TEST*/ 209