00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015 #include <parmetislib.h>
00016
00017
00018 #define DEBUG_IPART_
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028 void Moc_InitPartition_RB(CtrlType *ctrl, GraphType *graph, WorkSpaceType *wspace)
00029 {
00030 int i, j;
00031 int ncon, mype, npes, gnvtxs, ngroups;
00032 idxtype *xadj, *adjncy, *adjwgt, *vwgt;
00033 idxtype *part, *gwhere0, *gwhere1;
00034 idxtype *tmpwhere, *tmpvwgt, *tmpxadj, *tmpadjncy, *tmpadjwgt;
00035 GraphType *agraph;
00036 int lnparts, fpart, fpe, lnpes;
00037 int twoparts=2, numflag = 0, wgtflag = 3, moptions[10], edgecut, max_cut;
00038 floattype *mytpwgts, mytpwgts2[2], lbvec[MAXNCON], lbsum, min_lbsum, wsum;
00039 MPI_Comm ipcomm;
00040 struct {
00041 floattype sum;
00042 int rank;
00043 } lpesum, gpesum;
00044
00045 ncon = graph->ncon;
00046 ngroups = amax(amin(RIP_SPLIT_FACTOR, ctrl->npes), 1);
00047
00048 IFSET(ctrl->dbglvl, DBG_TIME, MPI_Barrier(ctrl->comm));
00049 IFSET(ctrl->dbglvl, DBG_TIME, starttimer(ctrl->InitPartTmr));
00050
00051 agraph = Moc_AssembleAdaptiveGraph(ctrl, graph, wspace);
00052 part = idxmalloc(agraph->nvtxs, "Moc_IP_RB: part");
00053 xadj = idxmalloc(agraph->nvtxs+1, "Moc_IP_RB: xadj");
00054 adjncy = idxmalloc(agraph->nedges, "Moc_IP_RB: adjncy");
00055 adjwgt = idxmalloc(agraph->nedges, "Moc_IP_RB: adjwgt");
00056 vwgt = idxmalloc(agraph->nvtxs*ncon, "Moc_IP_RB: vwgt");
00057
00058 idxcopy(agraph->nvtxs*ncon, agraph->vwgt, vwgt);
00059 idxcopy(agraph->nvtxs+1, agraph->xadj, xadj);
00060 idxcopy(agraph->nedges, agraph->adjncy, adjncy);
00061 idxcopy(agraph->nedges, agraph->adjwgt, adjwgt);
00062
00063 MPI_Comm_split(ctrl->gcomm, ctrl->mype % ngroups, 0, &ipcomm);
00064 MPI_Comm_rank(ipcomm, &mype);
00065 MPI_Comm_size(ipcomm, &npes);
00066
00067 gnvtxs = agraph->nvtxs;
00068
00069 gwhere0 = idxsmalloc(gnvtxs, 0, "Moc_IP_RB: gwhere0");
00070 gwhere1 = idxmalloc(gnvtxs, "Moc_IP_RB: gwhere1");
00071
00072
00073
00074 mytpwgts = fsmalloc(ctrl->nparts, 0.0, "mytpwgts");
00075 for (i=0; i<ctrl->nparts; i++)
00076 for (j=0; j<ncon; j++)
00077 mytpwgts[i] += ctrl->tpwgts[i*ncon+j];
00078 for (i=0; i<ctrl->nparts; i++)
00079 mytpwgts[i] /= (floattype)ncon;
00080
00081
00082
00083 moptions[0] = 0;
00084 moptions[7] = ctrl->sync + (ctrl->mype % ngroups) + 1;
00085
00086 lnparts = ctrl->nparts;
00087 fpart = fpe = 0;
00088 lnpes = npes;
00089 while (lnpes > 1 && lnparts > 1) {
00090
00091 mytpwgts2[0] = ssum(lnparts/2, mytpwgts+fpart);
00092 mytpwgts2[1] = 1.0-mytpwgts2[0];
00093
00094 if (ncon == 1)
00095 METIS_WPartGraphKway2(&agraph->nvtxs, agraph->xadj, agraph->adjncy,
00096 agraph->vwgt, agraph->adjwgt, &wgtflag, &numflag, &twoparts, mytpwgts2,
00097 moptions, &edgecut, part);
00098 else {
00099 METIS_mCPartGraphRecursive2(&agraph->nvtxs, &ncon, agraph->xadj,
00100 agraph->adjncy, agraph->vwgt, agraph->adjwgt, &wgtflag, &numflag,
00101 &twoparts, mytpwgts2, moptions, &edgecut, part);
00102 }
00103
00104 wsum = ssum(lnparts/2, mytpwgts+fpart);
00105 sscale(lnparts/2, 1.0/wsum, mytpwgts+fpart);
00106 sscale(lnparts-lnparts/2, 1.0/(1.0-wsum), mytpwgts+fpart+lnparts/2);
00107
00108
00109 if (mype < fpe+lnpes/2) {
00110 Moc_KeepPart(agraph, wspace, part, 0);
00111 lnpes = lnpes/2;
00112 lnparts = lnparts/2;
00113 }
00114 else {
00115 Moc_KeepPart(agraph, wspace, part, 1);
00116 fpart = fpart + lnparts/2;
00117 fpe = fpe + lnpes/2;
00118 lnpes = lnpes - lnpes/2;
00119 lnparts = lnparts - lnparts/2;
00120 }
00121 }
00122
00123
00124 if (lnparts == 1) {
00125
00126 if (mype == fpe) {
00127 for (i=0; i<agraph->nvtxs; i++)
00128 gwhere0[agraph->label[i]] = fpart;
00129 }
00130 }
00131
00132 else {
00133 if (ncon == 1)
00134 METIS_WPartGraphKway2(&agraph->nvtxs, agraph->xadj, agraph->adjncy,
00135 agraph->vwgt, agraph->adjwgt, &wgtflag, &numflag, &lnparts, mytpwgts+fpart,
00136 moptions, &edgecut, part);
00137 else
00138 METIS_mCPartGraphRecursive2(&agraph->nvtxs, &ncon, agraph->xadj,
00139 agraph->adjncy, agraph->vwgt, agraph->adjwgt, &wgtflag, &numflag,
00140 &lnparts, mytpwgts+fpart, moptions, &edgecut, part);
00141
00142 for (i=0; i<agraph->nvtxs; i++)
00143 gwhere0[agraph->label[i]] = fpart + part[i];
00144 }
00145
00146 MPI_Allreduce((void *)gwhere0, (void *)gwhere1, gnvtxs, IDX_DATATYPE, MPI_SUM, ipcomm);
00147
00148 if (ngroups > 1) {
00149 tmpxadj = agraph->xadj;
00150 tmpadjncy = agraph->adjncy;
00151 tmpadjwgt = agraph->adjwgt;
00152 tmpvwgt = agraph->vwgt;
00153 tmpwhere = agraph->where;
00154 agraph->xadj = xadj;
00155 agraph->adjncy = adjncy;
00156 agraph->adjwgt = adjwgt;
00157 agraph->vwgt = vwgt;
00158 agraph->where = gwhere1;
00159 agraph->vwgt = vwgt;
00160 agraph->nvtxs = gnvtxs;
00161 Moc_ComputeSerialBalance(ctrl, agraph, gwhere1, lbvec);
00162 lbsum = ssum(ncon, lbvec);
00163
00164 edgecut = ComputeSerialEdgeCut(agraph);
00165 MPI_Allreduce((void *)&edgecut, (void *)&max_cut, 1, MPI_INT, MPI_MAX, ctrl->gcomm);
00166 MPI_Allreduce((void *)&lbsum, (void *)&min_lbsum, 1, MPI_DOUBLE, MPI_MIN, ctrl->gcomm);
00167
00168 lpesum.sum = lbsum;
00169 if (min_lbsum < UNBALANCE_FRACTION * (floattype)(ncon)) {
00170 if (lbsum < UNBALANCE_FRACTION * (floattype)(ncon))
00171 lpesum.sum = (floattype) (edgecut);
00172 else
00173 lpesum.sum = (floattype) (max_cut);
00174 }
00175
00176 MPI_Comm_rank(ctrl->gcomm, &(lpesum.rank));
00177 MPI_Allreduce((void *)&lpesum, (void *)&gpesum, 1, MPI_DOUBLE_INT, MPI_MINLOC, ctrl->gcomm);
00178 MPI_Bcast((void *)gwhere1, gnvtxs, IDX_DATATYPE, gpesum.rank, ctrl->gcomm);
00179
00180 agraph->xadj = tmpxadj;
00181 agraph->adjncy = tmpadjncy;
00182 agraph->adjwgt = tmpadjwgt;
00183 agraph->vwgt = tmpvwgt;
00184 agraph->where = tmpwhere;
00185 }
00186
00187 idxcopy(graph->nvtxs, gwhere1+graph->vtxdist[ctrl->mype], graph->where);
00188
00189 FreeGraph(agraph);
00190 MPI_Comm_free(&ipcomm);
00191 GKfree((void **)&gwhere0, (void **)&gwhere1, (void **)&mytpwgts, (void **)&part, (void **)&xadj, (void **)&adjncy, (void **)&adjwgt, (void **)&vwgt, LTERM);
00192
00193 IFSET(ctrl->dbglvl, DBG_TIME, MPI_Barrier(ctrl->comm));
00194 IFSET(ctrl->dbglvl, DBG_TIME, stoptimer(ctrl->InitPartTmr));
00195
00196 }
00197
00198
00199
00200
00201
00202 void Moc_KeepPart(GraphType *graph, WorkSpaceType *wspace, idxtype *part, int mypart)
00203 {
00204 int h, i, j, k;
00205 int nvtxs, ncon, mynvtxs, mynedges;
00206 idxtype *xadj, *vwgt, *adjncy, *adjwgt, *label;
00207 idxtype *rename;
00208
00209 nvtxs = graph->nvtxs;
00210 ncon = graph->ncon;
00211 xadj = graph->xadj;
00212 vwgt = graph->vwgt;
00213 adjncy = graph->adjncy;
00214 adjwgt = graph->adjwgt;
00215 label = graph->label;
00216
00217 rename = idxmalloc(nvtxs, "Moc_KeepPart: rename");
00218
00219 for (mynvtxs=0, i=0; i<nvtxs; i++) {
00220 if (part[i] == mypart)
00221 rename[i] = mynvtxs++;
00222 }
00223
00224 for (mynvtxs=0, mynedges=0, j=xadj[0], i=0; i<nvtxs; i++) {
00225 if (part[i] == mypart) {
00226 for (; j<xadj[i+1]; j++) {
00227 k = adjncy[j];
00228 if (part[k] == mypart) {
00229 adjncy[mynedges] = rename[k];
00230 adjwgt[mynedges++] = adjwgt[j];
00231 }
00232 }
00233 j = xadj[i+1];
00234
00235 for (h=0; h<ncon; h++)
00236 vwgt[mynvtxs*ncon+h] = vwgt[i*ncon+h];
00237 label[mynvtxs] = label[i];
00238 xadj[++mynvtxs] = mynedges;
00239
00240 }
00241 else {
00242 j = xadj[i+1];
00243 }
00244 }
00245
00246 graph->nvtxs = mynvtxs;
00247 graph->nedges = mynedges;
00248
00249 free(rename);
00250 }
00251
00252