1 /*
2 * pspases.c
3 *
4 * This file contains ordering routines that are to be used with the
5 * parallel Cholesky factorization code PSPASES
6 *
7 * Started 10/14/97
8 * George
9 *
10 * $Id: pspases.c 10535 2011-07-11 04:29:44Z karypis $
11 *
12 */
13
14 #include <parmetislib.h>
15
16
17 /***********************************************************************************
18 * This function is the entry point of the serial ordering algorithm.
19 ************************************************************************************/
ParMETIS_SerialNodeND(idx_t * vtxdist,idx_t * xadj,idx_t * adjncy,idx_t * numflag,idx_t * options,idx_t * order,idx_t * sizes,MPI_Comm * comm)20 int ParMETIS_SerialNodeND(idx_t *vtxdist, idx_t *xadj, idx_t *adjncy,
21 idx_t *numflag, idx_t *options, idx_t *order, idx_t *sizes,
22 MPI_Comm *comm)
23 {
24 idx_t i, npes, mype;
25 ctrl_t *ctrl=NULL;
26 graph_t *agraph=NULL;
27 idx_t *perm=NULL, *iperm=NULL;
28 idx_t *sendcount, *displs;
29
30 /* Setup the ctrl */
31 ctrl = SetupCtrl(PARMETIS_OP_OMETIS, options, 1, 1, NULL, NULL, *comm);
32 npes = ctrl->npes;
33 mype = ctrl->mype;
34
35 if (!ispow2(npes)) {
36 if (mype == 0)
37 printf("Error: The number of processors must be a power of 2!\n");
38 FreeCtrl(&ctrl);
39 return METIS_ERROR;
40 }
41
42
43 if (*numflag > 0)
44 ChangeNumbering(vtxdist, xadj, adjncy, order, npes, mype, 1);
45
46 STARTTIMER(ctrl, ctrl->TotalTmr);
47 STARTTIMER(ctrl, ctrl->MoveTmr);
48
49 agraph = AssembleEntireGraph(ctrl, vtxdist, xadj, adjncy);
50
51 STOPTIMER(ctrl, ctrl->MoveTmr);
52
53 if (mype == 0) {
54 perm = imalloc(agraph->nvtxs, "PAROMETISS: perm");
55 iperm = imalloc(agraph->nvtxs, "PAROMETISS: iperm");
56
57 METIS_NodeNDP(agraph->nvtxs, agraph->xadj, agraph->adjncy,
58 agraph->vwgt, npes, NULL, perm, iperm, sizes);
59 }
60
61 STARTTIMER(ctrl, ctrl->MoveTmr);
62
63 /* Broadcast the sizes array */
64 gkMPI_Bcast((void *)sizes, 2*npes, IDX_T, 0, ctrl->gcomm);
65
66 /* Scatter the iperm */
67 sendcount = imalloc(npes, "PAROMETISS: sendcount");
68 displs = imalloc(npes, "PAROMETISS: displs");
69 for (i=0; i<npes; i++) {
70 sendcount[i] = vtxdist[i+1]-vtxdist[i];
71 displs[i] = vtxdist[i];
72 }
73
74 gkMPI_Scatterv((void *)iperm, sendcount, displs, IDX_T, (void *)order,
75 vtxdist[mype+1]-vtxdist[mype], IDX_T, 0, ctrl->gcomm);
76
77 STOPTIMER(ctrl, ctrl->MoveTmr);
78 STOPTIMER(ctrl, ctrl->TotalTmr);
79 IFSET(ctrl->dbglvl, DBG_TIME, PrintTimingInfo(ctrl));
80 IFSET(ctrl->dbglvl, DBG_TIME, gkMPI_Barrier(ctrl->gcomm));
81
82 gk_free((void **)&agraph->xadj, &agraph->adjncy, &perm, &iperm,
83 &sendcount, &displs, &agraph, LTERM);
84
85 if (*numflag > 0)
86 ChangeNumbering(vtxdist, xadj, adjncy, order, npes, mype, 0);
87
88 goto DONE;
89
90 DONE:
91 FreeCtrl(&ctrl);
92 return METIS_OK;
93 }
94
95
96
97 /*************************************************************************
98 * This function assembles the graph into a single processor
99 **************************************************************************/
AssembleEntireGraph(ctrl_t * ctrl,idx_t * vtxdist,idx_t * xadj,idx_t * adjncy)100 graph_t *AssembleEntireGraph(ctrl_t *ctrl, idx_t *vtxdist, idx_t *xadj, idx_t *adjncy)
101 {
102 idx_t i, gnvtxs, nvtxs, gnedges, nedges;
103 idx_t npes = ctrl->npes, mype = ctrl->mype;
104 idx_t *axadj, *aadjncy;
105 idx_t *recvcounts, *displs;
106 graph_t *agraph;
107
108 gnvtxs = vtxdist[npes];
109 nvtxs = vtxdist[mype+1]-vtxdist[mype];
110 nedges = xadj[nvtxs];
111
112 recvcounts = imalloc(npes, "AssembleGraph: recvcounts");
113 displs = imalloc(npes+1, "AssembleGraph: displs");
114
115 /* Gather all the xadj arrays first */
116 for (i=0; i<nvtxs; i++)
117 xadj[i] = xadj[i+1]-xadj[i];
118
119 axadj = imalloc(gnvtxs+1, "AssembleEntireGraph: axadj");
120
121 for (i=0; i<npes; i++) {
122 recvcounts[i] = vtxdist[i+1]-vtxdist[i];
123 displs[i] = vtxdist[i];
124 }
125
126 /* Assemble the xadj and then the adjncy */
127 gkMPI_Gatherv((void *)xadj, nvtxs, IDX_T, axadj, recvcounts, displs,
128 IDX_T, 0, ctrl->comm);
129
130 MAKECSR(i, nvtxs, xadj);
131 MAKECSR(i, gnvtxs, axadj);
132
133 /* Gather all the adjncy arrays next */
134 /* Determine the # of edges stored at each processor */
135 gkMPI_Allgather((void *)(&nedges), 1, IDX_T, (void *)recvcounts, 1, IDX_T, ctrl->comm);
136
137 displs[0] = 0;
138 for (i=1; i<npes+1; i++)
139 displs[i] = displs[i-1] + recvcounts[i-1];
140 gnedges = displs[npes];
141
142 aadjncy = imalloc(gnedges, "AssembleEntireGraph: aadjncy");
143
144 /* Assemble the xadj and then the adjncy */
145 gkMPI_Gatherv((void *)adjncy, nedges, IDX_T, aadjncy, recvcounts, displs, IDX_T, 0, ctrl->comm);
146
147 /* myprintf(ctrl, "Gnvtxs: %"PRIDX", Gnedges: %"PRIDX"\n", gnvtxs, gnedges); */
148
149 agraph = CreateGraph();
150 agraph->nvtxs = gnvtxs;
151 agraph->nedges = gnedges;
152 agraph->xadj = axadj;
153 agraph->adjncy = aadjncy;
154
155 return agraph;
156 }
157