1*abe9303eSLisandro Dalcin #include <petscmat.h> 2*abe9303eSLisandro Dalcin #include <petsc/private/partitionerimpl.h> /*I "petscpartitioner.h" I*/ 3*abe9303eSLisandro Dalcin 4*abe9303eSLisandro Dalcin typedef struct { 5*abe9303eSLisandro Dalcin MatPartitioning mp; 6*abe9303eSLisandro Dalcin } PetscPartitioner_MatPartitioning; 7*abe9303eSLisandro Dalcin 8*abe9303eSLisandro Dalcin static PetscErrorCode PetscPartitionerMatPartitioningGetMatPartitioning_MatPartitioning(PetscPartitioner part, MatPartitioning *mp) 9*abe9303eSLisandro Dalcin { 10*abe9303eSLisandro Dalcin PetscPartitioner_MatPartitioning *p = (PetscPartitioner_MatPartitioning *) part->data; 11*abe9303eSLisandro Dalcin 12*abe9303eSLisandro Dalcin PetscFunctionBegin; 13*abe9303eSLisandro Dalcin *mp = p->mp; 14*abe9303eSLisandro Dalcin PetscFunctionReturn(0); 15*abe9303eSLisandro Dalcin } 16*abe9303eSLisandro Dalcin 17*abe9303eSLisandro Dalcin /*@C 18*abe9303eSLisandro Dalcin PetscPartitionerMatPartitioningGetMatPartitioning - Get a MatPartitioning instance wrapped by this PetscPartitioner. 19*abe9303eSLisandro Dalcin 20*abe9303eSLisandro Dalcin Not Collective 21*abe9303eSLisandro Dalcin 22*abe9303eSLisandro Dalcin Input Parameters: 23*abe9303eSLisandro Dalcin . part - The PetscPartitioner 24*abe9303eSLisandro Dalcin 25*abe9303eSLisandro Dalcin Output Parameters: 26*abe9303eSLisandro Dalcin . mp - The MatPartitioning 27*abe9303eSLisandro Dalcin 28*abe9303eSLisandro Dalcin Level: developer 29*abe9303eSLisandro Dalcin 30*abe9303eSLisandro Dalcin .seealso DMPlexDistribute(), PetscPartitionerCreate() 31*abe9303eSLisandro Dalcin @*/ 32*abe9303eSLisandro Dalcin PetscErrorCode PetscPartitionerMatPartitioningGetMatPartitioning(PetscPartitioner part, MatPartitioning *mp) 33*abe9303eSLisandro Dalcin { 34*abe9303eSLisandro Dalcin PetscErrorCode ierr; 35*abe9303eSLisandro Dalcin 36*abe9303eSLisandro Dalcin PetscFunctionBegin; 37*abe9303eSLisandro Dalcin PetscValidHeaderSpecific(part, PETSCPARTITIONER_CLASSID, 1); 38*abe9303eSLisandro Dalcin PetscValidPointer(mp, 2); 39*abe9303eSLisandro Dalcin ierr = PetscUseMethod(part,"PetscPartitionerMatPartitioningGetMatPartitioning_C",(PetscPartitioner,MatPartitioning*),(part,mp));CHKERRQ(ierr); 40*abe9303eSLisandro Dalcin PetscFunctionReturn(0); 41*abe9303eSLisandro Dalcin } 42*abe9303eSLisandro Dalcin 43*abe9303eSLisandro Dalcin static PetscErrorCode PetscPartitionerDestroy_MatPartitioning(PetscPartitioner part) 44*abe9303eSLisandro Dalcin { 45*abe9303eSLisandro Dalcin PetscPartitioner_MatPartitioning *p = (PetscPartitioner_MatPartitioning *) part->data; 46*abe9303eSLisandro Dalcin PetscErrorCode ierr; 47*abe9303eSLisandro Dalcin 48*abe9303eSLisandro Dalcin PetscFunctionBegin; 49*abe9303eSLisandro Dalcin ierr = MatPartitioningDestroy(&p->mp);CHKERRQ(ierr); 50*abe9303eSLisandro Dalcin ierr = PetscFree(part->data);CHKERRQ(ierr); 51*abe9303eSLisandro Dalcin PetscFunctionReturn(0); 52*abe9303eSLisandro Dalcin } 53*abe9303eSLisandro Dalcin 54*abe9303eSLisandro Dalcin static PetscErrorCode PetscPartitionerView_MatPartitioning_ASCII(PetscPartitioner part, PetscViewer viewer) 55*abe9303eSLisandro Dalcin { 56*abe9303eSLisandro Dalcin PetscPartitioner_MatPartitioning *p = (PetscPartitioner_MatPartitioning *) part->data; 57*abe9303eSLisandro Dalcin PetscViewerFormat format; 58*abe9303eSLisandro Dalcin PetscErrorCode ierr; 59*abe9303eSLisandro Dalcin 60*abe9303eSLisandro Dalcin PetscFunctionBegin; 61*abe9303eSLisandro Dalcin ierr = PetscViewerGetFormat(viewer, &format);CHKERRQ(ierr); 62*abe9303eSLisandro Dalcin ierr = PetscViewerASCIIPrintf(viewer, "MatPartitioning Graph Partitioner:\n");CHKERRQ(ierr); 63*abe9303eSLisandro Dalcin ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 64*abe9303eSLisandro Dalcin if (p->mp) {ierr = MatPartitioningView(p->mp, viewer);CHKERRQ(ierr);} 65*abe9303eSLisandro Dalcin ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 66*abe9303eSLisandro Dalcin PetscFunctionReturn(0); 67*abe9303eSLisandro Dalcin } 68*abe9303eSLisandro Dalcin 69*abe9303eSLisandro Dalcin static PetscErrorCode PetscPartitionerView_MatPartitioning(PetscPartitioner part, PetscViewer viewer) 70*abe9303eSLisandro Dalcin { 71*abe9303eSLisandro Dalcin PetscBool iascii; 72*abe9303eSLisandro Dalcin PetscErrorCode ierr; 73*abe9303eSLisandro Dalcin 74*abe9303eSLisandro Dalcin PetscFunctionBegin; 75*abe9303eSLisandro Dalcin PetscValidHeaderSpecific(part, PETSCPARTITIONER_CLASSID, 1); 76*abe9303eSLisandro Dalcin PetscValidHeaderSpecific(viewer, PETSC_VIEWER_CLASSID, 2); 77*abe9303eSLisandro Dalcin ierr = PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);CHKERRQ(ierr); 78*abe9303eSLisandro Dalcin if (iascii) {ierr = PetscPartitionerView_MatPartitioning_ASCII(part, viewer);CHKERRQ(ierr);} 79*abe9303eSLisandro Dalcin PetscFunctionReturn(0); 80*abe9303eSLisandro Dalcin } 81*abe9303eSLisandro Dalcin 82*abe9303eSLisandro Dalcin static PetscErrorCode PetscPartitionerSetFromOptions_MatPartitioning(PetscOptionItems *PetscOptionsObject, PetscPartitioner part) 83*abe9303eSLisandro Dalcin { 84*abe9303eSLisandro Dalcin PetscPartitioner_MatPartitioning *p = (PetscPartitioner_MatPartitioning *) part->data; 85*abe9303eSLisandro Dalcin PetscErrorCode ierr; 86*abe9303eSLisandro Dalcin 87*abe9303eSLisandro Dalcin PetscFunctionBegin; 88*abe9303eSLisandro Dalcin ierr = PetscObjectSetOptionsPrefix((PetscObject)p->mp,((PetscObject)part)->prefix);CHKERRQ(ierr); 89*abe9303eSLisandro Dalcin ierr = MatPartitioningSetFromOptions(p->mp);CHKERRQ(ierr); 90*abe9303eSLisandro Dalcin PetscFunctionReturn(0); 91*abe9303eSLisandro Dalcin } 92*abe9303eSLisandro Dalcin 93*abe9303eSLisandro Dalcin static PetscErrorCode PetscPartitionerPartition_MatPartitioning(PetscPartitioner part, PetscInt nparts, PetscInt numVertices, PetscInt start[], PetscInt adjacency[], PetscSection vertSection, PetscSection targetSection, PetscSection partSection, IS *is) 94*abe9303eSLisandro Dalcin { 95*abe9303eSLisandro Dalcin PetscPartitioner_MatPartitioning *p = (PetscPartitioner_MatPartitioning *) part->data; 96*abe9303eSLisandro Dalcin Mat matadj; 97*abe9303eSLisandro Dalcin IS is1, is2, is3; 98*abe9303eSLisandro Dalcin PetscReal *tpwgts = NULL; 99*abe9303eSLisandro Dalcin PetscInt numVerticesGlobal, numEdges; 100*abe9303eSLisandro Dalcin PetscInt *i, *j, *vwgt = NULL; 101*abe9303eSLisandro Dalcin MPI_Comm comm; 102*abe9303eSLisandro Dalcin PetscErrorCode ierr; 103*abe9303eSLisandro Dalcin 104*abe9303eSLisandro Dalcin PetscFunctionBegin; 105*abe9303eSLisandro Dalcin ierr = PetscObjectGetComm((PetscObject)part, &comm);CHKERRQ(ierr); 106*abe9303eSLisandro Dalcin 107*abe9303eSLisandro Dalcin /* TODO: MatCreateMPIAdj should maybe take global number of ROWS */ 108*abe9303eSLisandro Dalcin /* TODO: And vertex distribution in PetscPartitionerPartition_ParMetis should be done using PetscSplitOwnership */ 109*abe9303eSLisandro Dalcin numVerticesGlobal = PETSC_DECIDE; 110*abe9303eSLisandro Dalcin ierr = PetscSplitOwnership(comm, &numVertices, &numVerticesGlobal);CHKERRQ(ierr); 111*abe9303eSLisandro Dalcin 112*abe9303eSLisandro Dalcin /* copy arrays to avoid memory errors because MatMPIAdjSetPreallocation copies just pointers */ 113*abe9303eSLisandro Dalcin numEdges = start[numVertices]; 114*abe9303eSLisandro Dalcin ierr = PetscMalloc1(numVertices+1, &i);CHKERRQ(ierr); 115*abe9303eSLisandro Dalcin ierr = PetscMalloc1(numEdges, &j);CHKERRQ(ierr); 116*abe9303eSLisandro Dalcin ierr = PetscArraycpy(i, start, numVertices+1);CHKERRQ(ierr); 117*abe9303eSLisandro Dalcin ierr = PetscArraycpy(j, adjacency, numEdges);CHKERRQ(ierr); 118*abe9303eSLisandro Dalcin 119*abe9303eSLisandro Dalcin /* construct the adjacency matrix */ 120*abe9303eSLisandro Dalcin ierr = MatCreateMPIAdj(comm, numVertices, numVerticesGlobal, i, j, NULL, &matadj);CHKERRQ(ierr); 121*abe9303eSLisandro Dalcin ierr = MatPartitioningSetAdjacency(p->mp, matadj);CHKERRQ(ierr); 122*abe9303eSLisandro Dalcin ierr = MatPartitioningSetNParts(p->mp, nparts);CHKERRQ(ierr); 123*abe9303eSLisandro Dalcin 124*abe9303eSLisandro Dalcin /* calculate partition weights */ 125*abe9303eSLisandro Dalcin if (targetSection) { 126*abe9303eSLisandro Dalcin PetscReal sumt; 127*abe9303eSLisandro Dalcin PetscInt p; 128*abe9303eSLisandro Dalcin 129*abe9303eSLisandro Dalcin sumt = 0.0; 130*abe9303eSLisandro Dalcin ierr = PetscMalloc1(nparts,&tpwgts);CHKERRQ(ierr); 131*abe9303eSLisandro Dalcin for (p = 0; p < nparts; ++p) { 132*abe9303eSLisandro Dalcin PetscInt tpd; 133*abe9303eSLisandro Dalcin 134*abe9303eSLisandro Dalcin ierr = PetscSectionGetDof(targetSection,p,&tpd);CHKERRQ(ierr); 135*abe9303eSLisandro Dalcin sumt += tpd; 136*abe9303eSLisandro Dalcin tpwgts[p] = tpd; 137*abe9303eSLisandro Dalcin } 138*abe9303eSLisandro Dalcin if (sumt) { /* METIS/ParMETIS do not like exactly zero weight */ 139*abe9303eSLisandro Dalcin for (p = 0, sumt = 0.0; p < nparts; ++p) { 140*abe9303eSLisandro Dalcin tpwgts[p] = PetscMax(tpwgts[p],PETSC_SMALL); 141*abe9303eSLisandro Dalcin sumt += tpwgts[p]; 142*abe9303eSLisandro Dalcin } 143*abe9303eSLisandro Dalcin for (p = 0; p < nparts; ++p) tpwgts[p] /= sumt; 144*abe9303eSLisandro Dalcin for (p = 0, sumt = 0.0; p < nparts-1; ++p) sumt += tpwgts[p]; 145*abe9303eSLisandro Dalcin tpwgts[nparts - 1] = 1. - sumt; 146*abe9303eSLisandro Dalcin } else { 147*abe9303eSLisandro Dalcin ierr = PetscFree(tpwgts);CHKERRQ(ierr); 148*abe9303eSLisandro Dalcin } 149*abe9303eSLisandro Dalcin } 150*abe9303eSLisandro Dalcin ierr = MatPartitioningSetPartitionWeights(p->mp, tpwgts);CHKERRQ(ierr); 151*abe9303eSLisandro Dalcin 152*abe9303eSLisandro Dalcin /* calculate vertex weights */ 153*abe9303eSLisandro Dalcin if (vertSection) { 154*abe9303eSLisandro Dalcin PetscInt v; 155*abe9303eSLisandro Dalcin 156*abe9303eSLisandro Dalcin ierr = PetscMalloc1(numVertices,&vwgt);CHKERRQ(ierr); 157*abe9303eSLisandro Dalcin for (v = 0; v < numVertices; ++v) { 158*abe9303eSLisandro Dalcin ierr = PetscSectionGetDof(vertSection, v, &vwgt[v]);CHKERRQ(ierr); 159*abe9303eSLisandro Dalcin } 160*abe9303eSLisandro Dalcin } 161*abe9303eSLisandro Dalcin ierr = MatPartitioningSetVertexWeights(p->mp, vwgt);CHKERRQ(ierr); 162*abe9303eSLisandro Dalcin 163*abe9303eSLisandro Dalcin /* apply the partitioning */ 164*abe9303eSLisandro Dalcin ierr = MatPartitioningApply(p->mp, &is1);CHKERRQ(ierr); 165*abe9303eSLisandro Dalcin 166*abe9303eSLisandro Dalcin /* construct the PetscSection */ 167*abe9303eSLisandro Dalcin { 168*abe9303eSLisandro Dalcin PetscInt v; 169*abe9303eSLisandro Dalcin const PetscInt *assignment_arr; 170*abe9303eSLisandro Dalcin 171*abe9303eSLisandro Dalcin ierr = ISGetIndices(is1, &assignment_arr);CHKERRQ(ierr); 172*abe9303eSLisandro Dalcin for (v = 0; v < numVertices; ++v) {ierr = PetscSectionAddDof(partSection, assignment_arr[v], 1);CHKERRQ(ierr);} 173*abe9303eSLisandro Dalcin ierr = ISRestoreIndices(is1, &assignment_arr);CHKERRQ(ierr); 174*abe9303eSLisandro Dalcin } 175*abe9303eSLisandro Dalcin 176*abe9303eSLisandro Dalcin /* convert assignment IS to global numbering IS */ 177*abe9303eSLisandro Dalcin ierr = ISPartitioningToNumbering(is1, &is2);CHKERRQ(ierr); 178*abe9303eSLisandro Dalcin ierr = ISDestroy(&is1);CHKERRQ(ierr); 179*abe9303eSLisandro Dalcin 180*abe9303eSLisandro Dalcin /* renumber IS into local numbering */ 181*abe9303eSLisandro Dalcin ierr = ISOnComm(is2, PETSC_COMM_SELF, PETSC_USE_POINTER, &is1);CHKERRQ(ierr); 182*abe9303eSLisandro Dalcin ierr = ISRenumber(is1, NULL, NULL, &is3);CHKERRQ(ierr); 183*abe9303eSLisandro Dalcin ierr = ISDestroy(&is1);CHKERRQ(ierr); 184*abe9303eSLisandro Dalcin ierr = ISDestroy(&is2);CHKERRQ(ierr); 185*abe9303eSLisandro Dalcin 186*abe9303eSLisandro Dalcin /* invert IS */ 187*abe9303eSLisandro Dalcin ierr = ISSetPermutation(is3);CHKERRQ(ierr); 188*abe9303eSLisandro Dalcin ierr = ISInvertPermutation(is3, numVertices, &is1);CHKERRQ(ierr); 189*abe9303eSLisandro Dalcin ierr = ISDestroy(&is3);CHKERRQ(ierr); 190*abe9303eSLisandro Dalcin 191*abe9303eSLisandro Dalcin ierr = MatDestroy(&matadj);CHKERRQ(ierr); 192*abe9303eSLisandro Dalcin *is = is1; 193*abe9303eSLisandro Dalcin PetscFunctionReturn(0); 194*abe9303eSLisandro Dalcin } 195*abe9303eSLisandro Dalcin 196*abe9303eSLisandro Dalcin static PetscErrorCode PetscPartitionerInitialize_MatPartitioning(PetscPartitioner part) 197*abe9303eSLisandro Dalcin { 198*abe9303eSLisandro Dalcin PetscErrorCode ierr; 199*abe9303eSLisandro Dalcin 200*abe9303eSLisandro Dalcin PetscFunctionBegin; 201*abe9303eSLisandro Dalcin part->ops->view = PetscPartitionerView_MatPartitioning; 202*abe9303eSLisandro Dalcin part->ops->setfromoptions = PetscPartitionerSetFromOptions_MatPartitioning; 203*abe9303eSLisandro Dalcin part->ops->destroy = PetscPartitionerDestroy_MatPartitioning; 204*abe9303eSLisandro Dalcin part->ops->partition = PetscPartitionerPartition_MatPartitioning; 205*abe9303eSLisandro Dalcin ierr = PetscObjectComposeFunction((PetscObject)part,"PetscPartitionerMatPartitioningGetMatPartitioning_C",PetscPartitionerMatPartitioningGetMatPartitioning_MatPartitioning);CHKERRQ(ierr); 206*abe9303eSLisandro Dalcin PetscFunctionReturn(0); 207*abe9303eSLisandro Dalcin } 208*abe9303eSLisandro Dalcin 209*abe9303eSLisandro Dalcin /*MC 210*abe9303eSLisandro Dalcin PETSCPARTITIONERMATPARTITIONING = "matpartitioning" - A PetscPartitioner object 211*abe9303eSLisandro Dalcin 212*abe9303eSLisandro Dalcin Level: developer 213*abe9303eSLisandro Dalcin 214*abe9303eSLisandro Dalcin .seealso: PetscPartitionerType, PetscPartitionerCreate(), PetscPartitionerSetType() 215*abe9303eSLisandro Dalcin M*/ 216*abe9303eSLisandro Dalcin 217*abe9303eSLisandro Dalcin PETSC_EXTERN PetscErrorCode PetscPartitionerCreate_MatPartitioning(PetscPartitioner part) 218*abe9303eSLisandro Dalcin { 219*abe9303eSLisandro Dalcin PetscPartitioner_MatPartitioning *p; 220*abe9303eSLisandro Dalcin PetscErrorCode ierr; 221*abe9303eSLisandro Dalcin 222*abe9303eSLisandro Dalcin PetscFunctionBegin; 223*abe9303eSLisandro Dalcin PetscValidHeaderSpecific(part, PETSCPARTITIONER_CLASSID, 1); 224*abe9303eSLisandro Dalcin ierr = PetscNewLog(part, &p);CHKERRQ(ierr); 225*abe9303eSLisandro Dalcin part->data = p; 226*abe9303eSLisandro Dalcin ierr = PetscPartitionerInitialize_MatPartitioning(part);CHKERRQ(ierr); 227*abe9303eSLisandro Dalcin ierr = MatPartitioningCreate(PetscObjectComm((PetscObject)part), &p->mp);CHKERRQ(ierr); 228*abe9303eSLisandro Dalcin PetscFunctionReturn(0); 229*abe9303eSLisandro Dalcin } 230*abe9303eSLisandro Dalcin 231