1 /* This software was developed by Bruce Hendrickson and Robert Leland   *
2  * at Sandia National Laboratories under US Department of Energy        *
3  * contract DE-AC04-76DP00789 and is copyrighted by Sandia Corporation. */
4 
5 #include <stdio.h>
6 #include "params.h"
7 #include "defs.h"
8 #include "structs.h"
9 
10 
11 int       DEBUG_COVER = 0;	/* Debug output in min vtx cover routines? */
12 int       FLATTEN = FALSE;	/* Merge indistinguishable vtxs first? */
13 int       COARSE_KLV = TRUE;	/* Use KLV as multilevel refinement? */
14 int       COARSE_BPM = FALSE;	/* Use vertex cover as ML refinement? */
15 int       VERTEX_COVER = 0;	/* Do vertex cover to refine separator? */
16 
17  /* Once, or iteratively? */
18 
coarsen_klv(graph,nvtxs,nedges,using_vwgts,using_ewgts,term_wgts,igeom,coords,vwgt_max,assignment,goal,architecture,hops,solver_flag,ndims,nsets,vmax,mediantype,mkconnected,eigtol,nstep,step,pbndy_list,weights,give_up)19 void      coarsen_klv(graph, nvtxs, nedges, using_vwgts, using_ewgts, term_wgts,
20 		      igeom, coords, vwgt_max, assignment, goal, architecture, hops,
21 		          solver_flag, ndims, nsets, vmax, mediantype, mkconnected,
22 		                eigtol, nstep, step, pbndy_list, weights, give_up)
23 /* Coarsen until nvtxs < vmax, compute and uncoarsen. */
24 struct vtx_data **graph;	/* array of vtx data for graph */
25 int       nvtxs;		/* number of vertices in graph */
26 int       nedges;		/* number of edges in graph */
27 int       using_vwgts;		/* are vertices weights being used? */
28 int       using_ewgts;		/* are edge weights being used? */
29 float    *term_wgts[];		/* weights for terminal propogation */
30 int       igeom;                /* dimension for geometric information */
31 float   **coords;               /* coordinates for vertices */
32 int       vwgt_max;		/* largest vertex weight */
33 short    *assignment;		/* processor each vertex gets assigned to */
34 double   *goal;			/* desired set sizes */
35 int       architecture;		/* 0 => hypercube, d => d-dimensional mesh */
36 short     (*hops)[MAXSETS];	/* cost of edge between sets */
37 int       solver_flag;		/* which eigensolver to use */
38 int       ndims;		/* number of eigenvectors to calculate */
39 int       nsets;		/* number of sets being divided into */
40 int       vmax;			/* largest subgraph to stop coarsening */
41 int       mediantype;		/* flag for different assignment strategies */
42 int       mkconnected;		/* enforce connectivity before eigensolve? */
43 double    eigtol;		/* tolerence in eigen calculation */
44 int       nstep;		/* number of coarsenings between RQI steps */
45 int       step;			/* current step number */
46 int     **pbndy_list;		/* pointer to returned boundary list */
47 double   *weights;		/* weights of vertices in each set */
48 int       give_up;		/* has coarsening bogged down? */
49 {
50     extern FILE *Output_File;	/* output file or null */
51     extern int DEBUG_TRACE;	/* trace the execution of the code */
52     extern int DEBUG_COARSEN;	/* debug flag for coarsening */
53     extern int DEBUG_CONNECTED;	/* debug flag for connectivity checking */
54     extern double COARSEN_RATIO_MIN;	/* vtx reduction demanded */
55     extern int COARSEN_VWGTS;	/* use vertex weights while coarsening? */
56     extern int COARSEN_EWGTS;	/* use edge weights while coarsening? */
57     extern int COARSE_KL_BOTTOM;/* force KL invocation at bottom level? */
58     extern int LIMIT_KL_EWGTS;	/* limit edges weights in KL? */
59     extern int COARSE_KLV;	/* apply klv as a smoother? */
60     extern int COARSE_BPM;	/* apply bipartite matching/flow as a smoother? */
61     extern double KL_IMBALANCE;	/* fractional imbalance allowed in KL */
62     struct connect_data *cdata;	/* data structure for enforcing connectivity */
63     struct vtx_data **cgraph;	/* array of vtx data for coarsened graph */
64     double   *yvecs[MAXDIMS + 1];	/* eigenvectors for subgraph */
65     double    evals[MAXDIMS + 1];	/* eigenvalues returned */
66     double    new_goal[MAXSETS];/* new goal if not using vertex weights */
67     double   *real_goal;	/* chooses between goal and new_goal */
68     double    goal_weight;	/* total weight of vertices in goal */
69     double   *vwsqrt = NULL;	/* square root of vertex weights */
70     double    maxdeg;		/* maximum weighted degree of a vertex */
71     float    *cterm_wgts[MAXSETS];	/* terminal weights for coarse graph */
72     float    *new_term_wgts[MAXSETS];	/* modified for Bui's method */
73     float   **real_term_wgts;	/* which of previous two to use */
74     float     ewgt_max;		/* largest edge weight in graph */
75     float    *twptr;		/* loops through term_wgts */
76     float    *twptr_save;	/* copy of twptr */
77     float    *ctwptr;		/* loops through cterm_wgts */
78     float   **ccoords;		/* coarse graph coordinates */
79     int      *active;		/* space for assign routine */
80     int      *v2cv;		/* mapping from vtxs to coarse vtxs */
81     int      *flag;		/* scatter array for coarse bndy vtxs */
82     int      *bndy_list;	/* list of vertices on boundary */
83     int      *cbndy_list;	/* list of vertices of coarse graph on boundary */
84     short    *cassignment;	/* set assignments for coarsened vertices */
85     int       flattened;	/* was this graph flattened? */
86     int       list_length;	/* length of boundary vtx list */
87     int       cnvtxs;		/* number of vertices in coarsened graph */
88     int       cnedges;		/* number of edges in coarsened graph */
89     int       cvwgt_max;	/* largest vertex weight in coarsened graph */
90     int       nextstep;		/* next step in RQI test */
91     int       max_dev;		/* largest allowed deviation from balance */
92     int       i, j;		/* loop counters */
93     int       sfree(), find_bndy(), flatten();
94     double   *smalloc(), find_maxdeg();
95     void      makevwsqrt(), make_connected(), print_connected(), eigensolve();
96     void      make_unconnected(), assign(), klvspiff(), coarsen1(), free_graph();
97     void      compress_ewgts(), restore_ewgts(), count_weights(), bpm_improve();
98 
99     if (DEBUG_COARSEN > 0 || DEBUG_TRACE > 0) {
100 	printf("<Entering coarsen_kl, step=%d, nvtxs=%d, nedges=%d, vmax=%d>\n",
101 	       step, nvtxs, nedges, vmax);
102     }
103 
104     /* Is problem small enough to solve? */
105     if (nvtxs <= vmax || give_up) {
106 
107 	/* Starting here, everything could be replaced by simple partitioner */
108 	/* coded w/ few lines below. */
109 /*
110 	if (using_vwgts) {
111 	    vwsqrt = (double *) smalloc((unsigned) (nvtxs + 1) * sizeof(double));
112 	    makevwsqrt(vwsqrt, graph, nvtxs);
113 	}
114 	else
115 	    vwsqrt = NULL;
116 
117 	/* Create space for subgraph yvecs. */
118 /*
119 	for (i = 1; i <= ndims; i++) {
120 	    yvecs[i] = (double *) smalloc((unsigned) (nvtxs + 1) * sizeof(double));
121 	}
122 
123 	active = (int *) smalloc((unsigned) nvtxs * sizeof(int));
124 	if (mkconnected) {
125 	    make_connected(graph, nvtxs, &nedges, (short *) &(yvecs[1][0]),
126 			   active, &cdata, using_ewgts);
127 	    if (DEBUG_CONNECTED > 0) {
128 		printf("Enforcing connectivity on coarse graph\n");
129 		print_connected(cdata);
130 	    }
131 	}
132 
133 	/* If not coarsening ewgts, then need care with term_wgts. */
134 /*
135 	if (!using_ewgts && term_wgts[1] != NULL && step != 0) {
136 	    twptr = (float *) smalloc((unsigned)
137 				      (nvtxs + 1) * (nsets - 1) * sizeof(float));
138 	    twptr_save = twptr;
139 	    for (j = 1; j < nsets; j++) {
140 		new_term_wgts[j] = twptr;
141 		twptr += nvtxs + 1;
142 	    }
143 
144 	    for (j = 1; j < nsets; j++) {
145 		twptr = term_wgts[j];
146 		ctwptr = new_term_wgts[j];
147 		for (i = 1; i <= nvtxs; i++) {
148 		    if (twptr[i] > .5)
149 			ctwptr[i] = 1;
150 		    else if (twptr[i] < -.5)
151 			ctwptr[i] = -1;
152 		    else
153 			ctwptr[i] = 0;
154 		}
155 	    }
156 	    real_term_wgts = new_term_wgts;
157 	}
158 	else {
159 	    real_term_wgts = term_wgts;
160 	    new_term_wgts[1] = NULL;
161 	}
162 
163 
164 	/* Partition coarse graph with spectral method */
165 /*
166 	maxdeg = find_maxdeg(graph, nvtxs, using_ewgts, &ewgt_max);
167 	eigensolve(graph, nvtxs, nedges, maxdeg, vwgt_max, vwsqrt,
168 		   using_vwgts, using_ewgts, real_term_wgts, 0, (float **) NULL,
169 		   yvecs, evals, architecture, assignment, goal,
170 		   solver_flag, FALSE, 0, ndims, mediantype, eigtol);
171 
172 	if (mkconnected)
173 	    make_unconnected(graph, &nedges, &cdata, using_ewgts);
174 
175 	if (!COARSEN_VWGTS && step != 0) {	/* Construct new goal */
176 /*
177 	    goal_weight = 0;
178 	    for (i = 0; i < nsets; i++)
179 		goal_weight += goal[i];
180 	    for (i = 0; i < nsets; i++)
181 		new_goal[i] = goal[i] * (nvtxs / goal_weight);
182 	    real_goal = new_goal;
183 	}
184 	else
185 	    real_goal = goal;
186 
187 	assign(graph, yvecs, nvtxs, ndims, architecture, nsets, vwsqrt, assignment,
188 	       active, mediantype, real_goal, vwgt_max);
189 
190 	sfree((char *) active);
191 
192 	/* Check for boundary and flag those guys with a 2 */
193 /*
194 	list_length = find_bndy(graph, nvtxs, assignment, 2, &bndy_list);
195 	count_weights(graph, nvtxs, assignment, nsets + 1, weights, (vwgt_max != 1));
196 
197 	if (COARSE_KL_BOTTOM || !(step % nstep)) {
198 	    if (LIMIT_KL_EWGTS) {
199 		find_maxdeg(graph, nvtxs, using_ewgts, &ewgt_max);
200 		compress_ewgts(graph, nvtxs, nedges, ewgt_max, using_ewgts);
201 	    }
202 
203 	    max_dev = (step == 0) ? vwgt_max : 5 * vwgt_max;
204 	    goal_weight = 0;
205 	    for (i = 0; i < nsets; i++) {
206 		goal_weight += real_goal[i];
207 	    }
208 	    goal_weight *= KL_IMBALANCE / nsets;
209 	    if (goal_weight > max_dev) {
210 		max_dev = goal_weight;
211 	    }
212 
213 	    if (COARSE_KLV) {
214 		klvspiff(graph, nvtxs, assignment, real_goal,
215 			 max_dev, &bndy_list, weights);
216 	    }
217 	    if (COARSE_BPM) {
218 		bpm_improve(graph, assignment, real_goal, max_dev, &bndy_list,
219 			    weights, using_vwgts);
220 	    }
221 	    if (LIMIT_KL_EWGTS)
222 		restore_ewgts(graph, nvtxs);
223 	}
224 
225 	if (real_term_wgts != term_wgts && new_term_wgts[1] != NULL) {
226 	    sfree((char *) real_term_wgts[1]);
227 	}
228 	if (vwsqrt != NULL)
229 	    sfree((char *) vwsqrt);
230 	for (i = ndims; i > 0; i--)
231 	    sfree((char *) yvecs[i]);
232 
233 	real_goal = goal;
234 */
235 
236 	/* The following could replace everything above. */
237 /* */
238 	real_goal = goal;
239 
240         simple_part(graph, nvtxs, assignment, nsets, 1, real_goal);
241         list_length = find_bndy(graph, nvtxs, assignment, 2, &bndy_list);
242 
243         count_weights(graph, nvtxs, assignment, nsets + 1, weights, 1);
244 
245         max_dev = (step == 0) ? vwgt_max : 5 * vwgt_max;
246         goal_weight = 0;
247         for (i = 0; i < nsets; i++)
248             goal_weight += real_goal[i];
249         goal_weight *= KL_IMBALANCE / nsets;
250         if (goal_weight > max_dev)
251             max_dev = goal_weight;
252 
253 	if (COARSE_KLV) {
254             klvspiff(graph, nvtxs, assignment, real_goal,
255                  max_dev, &bndy_list, weights);
256 	}
257 	if (COARSE_BPM) {
258 	    bpm_improve(graph, assignment, real_goal, max_dev, &bndy_list,
259 			    weights, using_vwgts);
260 	}
261 /* */
262 
263 	*pbndy_list = bndy_list;
264 	return;
265     }
266 
267     /* Otherwise I have to coarsen. */
268     flattened = FALSE;
269     if (coords != NULL) {
270         ccoords = (float **) smalloc((unsigned) igeom * sizeof(float *));
271     }
272     else {
273         ccoords = NULL;
274     }
275     if (FLATTEN && step == 0) {
276 	flattened = flatten(graph, nvtxs, nedges, &cgraph, &cnvtxs, &cnedges, &v2cv,
277 			    using_ewgts && COARSEN_EWGTS,
278 			    igeom, coords, ccoords);
279     }
280     if (!flattened) {
281 	coarsen1(graph, nvtxs, nedges, &cgraph, &cnvtxs, &cnedges, &v2cv,
282 		 igeom, coords, ccoords,
283 		 using_ewgts && COARSEN_EWGTS);
284     }
285 
286     if (term_wgts[1] != NULL) {
287 	twptr = (float *) smalloc((unsigned)
288 				  (cnvtxs + 1) * (nsets - 1) * sizeof(float));
289 	twptr_save = twptr;
290 	for (i = (cnvtxs + 1) * (nsets - 1); i; i--) {
291 	    *twptr++ = 0;
292 	}
293 	twptr = twptr_save;
294 	for (j = 1; j < nsets; j++) {
295 	    cterm_wgts[j] = twptr;
296 	    twptr += cnvtxs + 1;
297 	}
298 	for (j = 1; j < nsets; j++) {
299 	    ctwptr = cterm_wgts[j];
300 	    twptr = term_wgts[j];
301 	    for (i = 1; i < nvtxs; i++) {
302 		ctwptr[v2cv[i]] += twptr[i];
303 	    }
304 	}
305     }
306 
307     else {
308 	cterm_wgts[1] = NULL;
309     }
310 
311     /* If coarsening isn't working very well, give up and partition. */
312     give_up = FALSE;
313     if (nvtxs * COARSEN_RATIO_MIN < cnvtxs && cnvtxs > vmax && !flattened) {
314 	printf("WARNING: Coarsening not making enough progress, nvtxs = %d, cnvtxs = %d.\n",
315 	       nvtxs, cnvtxs);
316 	printf("         Recursive coarsening being stopped prematurely.\n");
317 	if (Output_File != NULL) {
318 	    fprintf(Output_File,
319 		    "WARNING: Coarsening not making enough progress, nvtxs = %d, cnvtxs = %d.\n",
320 		    nvtxs, cnvtxs);
321 	    fprintf(Output_File,
322 		    "         Recursive coarsening being stopped prematurely.\n");
323 	}
324 	give_up = TRUE;
325     }
326 
327     /* Now recurse on coarse subgraph. */
328     if (COARSEN_VWGTS) {
329 	cvwgt_max = 0;
330 	for (i = 1; i <= cnvtxs; i++) {
331 	    if (cgraph[i]->vwgt > cvwgt_max)
332 		cvwgt_max = cgraph[i]->vwgt;
333 	}
334     }
335 
336     else
337 	cvwgt_max = 1;
338 
339     cassignment = (short *) smalloc((unsigned) (cnvtxs + 1) * sizeof(short));
340     if (flattened)
341 	nextstep = step;
342     else
343 	nextstep = step + 1;
344     coarsen_klv(cgraph, cnvtxs, cnedges, COARSEN_VWGTS, COARSEN_EWGTS, cterm_wgts,
345 		igeom, ccoords, cvwgt_max, cassignment, goal, architecture, hops,
346 		solver_flag, ndims, nsets, vmax, mediantype, mkconnected,
347 		eigtol, nstep, nextstep, &cbndy_list, weights, give_up);
348 
349     /* Interpolate assignment back to fine graph. */
350     for (i = 1; i <= nvtxs; i++) {
351 	assignment[i] = cassignment[v2cv[i]];
352     }
353 
354     /* Construct boundary list from coarse boundary list. */
355     flag = (int *) smalloc((unsigned) (cnvtxs + 1) * sizeof(int));
356     for (i = 1; i <= cnvtxs; i++) {
357 	flag[i] = FALSE;
358     }
359     for (i = 0; cbndy_list[i]; i++) {
360 	flag[cbndy_list[i]] = TRUE;
361     }
362 
363     list_length = 0;
364     for (i = 1; i <= nvtxs; i++) {
365 	if (flag[v2cv[i]])
366 	    ++list_length;
367     }
368 
369     bndy_list = (int *) smalloc((unsigned) (list_length + 1) * sizeof(int));
370 
371     list_length = 0;
372     for (i = 1; i <= nvtxs; i++) {
373 	if (flag[v2cv[i]])
374 	    bndy_list[list_length++] = i;
375     }
376     bndy_list[list_length] = 0;
377 
378     sfree((char *) flag);
379     sfree((char *) cbndy_list);
380 
381 
382     /* Free the space that was allocated. */
383     sfree((char *) cassignment);
384     if (cterm_wgts[1] != NULL)
385 	sfree((char *) cterm_wgts[1]);
386     free_graph(cgraph);
387     sfree((char *) v2cv);
388 
389     /* Smooth using KL or BPM every nstep steps. */
390     if (!(step % nstep) && !flattened) {
391 /*
392     if (!(step % nstep)) {
393 */
394 	if (!COARSEN_VWGTS && step != 0) {	/* Construct new goal */
395 	    goal_weight = 0;
396 	    for (i = 0; i < nsets; i++)
397 		goal_weight += goal[i];
398 	    for (i = 0; i < nsets; i++)
399 		new_goal[i] = goal[i] * (nvtxs / goal_weight);
400 	    real_goal = new_goal;
401 	}
402 	else
403 	    real_goal = goal;
404 
405 	if (LIMIT_KL_EWGTS) {
406 	    find_maxdeg(graph, nvtxs, using_ewgts, &ewgt_max);
407 	    compress_ewgts(graph, nvtxs, nedges, ewgt_max,
408 			   using_ewgts);
409 	}
410 
411 	/* If not coarsening ewgts, then need care with term_wgts. */
412 	if (!using_ewgts && term_wgts[1] != NULL && step != 0) {
413 	    twptr = (float *) smalloc((unsigned)
414 				      (nvtxs + 1) * (nsets - 1) * sizeof(float));
415 	    twptr_save = twptr;
416 	    for (j = 1; j < nsets; j++) {
417 		new_term_wgts[j] = twptr;
418 		twptr += nvtxs + 1;
419 	    }
420 
421 	    for (j = 1; j < nsets; j++) {
422 		twptr = term_wgts[j];
423 		ctwptr = new_term_wgts[j];
424 		for (i = 1; i <= nvtxs; i++) {
425 		    if (twptr[i] > .5)
426 			ctwptr[i] = 1;
427 		    else if (twptr[i] < -.5)
428 			ctwptr[i] = -1;
429 		    else
430 			ctwptr[i] = 0;
431 		}
432 	    }
433 	    real_term_wgts = new_term_wgts;
434 	}
435 	else {
436 	    real_term_wgts = term_wgts;
437 	    new_term_wgts[1] = NULL;
438 	}
439 
440 	max_dev = (step == 0) ? vwgt_max : 5 * vwgt_max;
441 	goal_weight = 0;
442 	for (i = 0; i < nsets; i++) {
443 	    goal_weight += real_goal[i];
444 	}
445 	goal_weight *= KL_IMBALANCE / nsets;
446 	if (goal_weight > max_dev) {
447 	    max_dev = goal_weight;
448 	}
449 
450 	if (!COARSEN_VWGTS) {
451 	    count_weights(graph, nvtxs, assignment, nsets + 1, weights,
452 			  (vwgt_max != 1));
453 	}
454 
455 	if (COARSE_KLV) {
456 /*
457 printf("Before KLV");
458 print_sep_size(bndy_list, graph);
459 */
460 	    klvspiff(graph, nvtxs, assignment, real_goal,
461 		     max_dev, &bndy_list, weights);
462 	}
463 	if (COARSE_BPM) {
464 /*
465 printf("Before BPM");
466 print_sep_size(bndy_list, graph);
467 */
468 	    bpm_improve(graph, assignment, real_goal, max_dev, &bndy_list,
469 			weights, using_vwgts);
470 	}
471 /*
472 printf("Returning");
473 print_sep_size(bndy_list, graph);
474 */
475 
476 	if (real_term_wgts != term_wgts && new_term_wgts[1] != NULL) {
477 	    sfree((char *) real_term_wgts[1]);
478 	}
479 
480 	if (LIMIT_KL_EWGTS)
481 	    restore_ewgts(graph, nvtxs);
482     }
483 
484     *pbndy_list = bndy_list;
485 
486     /* Free the space that was allocated. */
487     if (ccoords != NULL) {
488         for (i = 0; i < igeom; i++)
489             sfree((char *) ccoords[i]);
490         sfree((char *) ccoords);
491     }
492 
493     if (DEBUG_COARSEN > 0) {
494 	printf(" Leaving coarsen_klv, step=%d\n", step);
495     }
496 }
497 
498 
print_sep_size(list,graph)499 void      print_sep_size(list, graph)
500 int      *list;
501 struct vtx_data **graph;	/* array of vtx data for graph */
502 {
503     int       sep_size, sep_weight;
504     int       i;
505 
506     sep_size = sep_weight = 0;
507     for (i = 0; list[i] != 0; i++) {
508 	sep_size++;
509 	sep_weight += graph[list[i]]->vwgt;
510     }
511     printf(" Sep_size = %d, Sep_weight = %d\n", sep_size, sep_weight);
512 }
513