1/*
2 * metis_util.hpp
3 *
4 * Created on: Nov 21, 2014
5 * Author: Pietro Incardona
6 */
7
8#ifndef METIS_UTIL_HPP
9#define METIS_UTIL_HPP
10
11#include <iostream>
12#include "metis.h"
13#include "SubdomainGraphNodes.hpp"
14#include "VTKWriter/VTKWriter.hpp"
15
16/*! \brief Metis graph structure
17 *
18 * Metis graph structure
19 *
20 */
21struct Metis_graph
22{
23 //! The number of vertices in the graph
24 idx_t * nvtxs;
25
26 //! number of balancing constrains
27 //! more practical, are the number of weights for each vertex
28 //! PS even we you specify vwgt == NULL ncon must be set at leat to
29 //! one
30 idx_t * ncon;
31
32 //! For each vertex it store the adjacency lost start for the vertex i
33 idx_t * xadj;
34
35 //! For each vertex it store a list of all neighborhood vertex
36 idx_t * adjncy;
37
38 //! Array that store the weight for each vertex
39 idx_t * vwgt;
40
41 //! Array of the vertex size, basically is the total communication amount
42 idx_t * vsize;
43 //! The weight of the edge
44 idx_t * adjwgt;
45 //! number of part to partition the graph
46 idx_t * nparts;
47 //! Desired weight for each partition (one for each constrain)
48 real_t * tpwgts;
49 //! For each partition load imbalance tollerated
50 real_t * ubvec;
51 //! Additional option for the graph partitioning
52 idx_t * options;
53 //! return the total comunication cost for each partition
54 idx_t * objval;
55 //! Is a output vector containing the partition for each vertex
56 idx_t * part;
57};
58
59//! Balance communication and computation
60#define BALANCE_CC 1
61//! Balance communication computation and memory
62#define BALANCE_CCM 2
63//! Balance computation and comunication and others
64#define BALANCE_CC_O(c) c+1
65
66/*! \brief Helper class to define Metis graph
67 *
68 * \tparam graph structure that store the graph
69 *
70 */
71
72template<typename Graph>
73class Metis
74{
75 //! Graph in metis reppresentation
76 Metis_graph Mg;
77
78 //! Original graph
79 Graph & g;
80
81 //! indicate how many time decompose/refine/re-decompose has been called
82 size_t n_dec;
83
84 //! Distribution tolerance
85 real_t dist_tol = 1.05;
86
87 /*! \brief Construct Adjacency list
88 *
89 * \param g Graph
90 *
91 */
92 void constructAdjList(Graph & g)
93 {
94 // create xadj and adjlist
95
96 // create xadj, adjlist, vwgt, adjwgt and vsize
97 Mg.xadj = new idx_t[g.getNVertex() + 1];
98 Mg.adjncy = new idx_t[g.getNEdge()];
99
100 //! starting point in the adjacency list
101 size_t prev = 0;
102
103 // actual position
104 size_t id = 0;
105
106 // for each vertex calculate the position of the starting point in the adjacency list
107 for (size_t i = 0; i < g.getNVertex(); i++)
108 {
109 // Calculate the starting point in the adjacency list
110 Mg.xadj[id] = prev;
111
112 // Create the adjacency list
113 for (size_t s = 0; s < g.getNChilds(i); s++)
114 {
115 Mg.adjncy[prev + s] = g.getChild(i, s);
116 }
117
118 // update the position for the next vertex
119 prev += g.getNChilds(i);
120
121 id++;
122 }
123
124 // Fill the last point
125 Mg.xadj[id] = prev;
126 }
127
128 /*! \brief Construct Adjacency list
129 *
130 * \param g Graph
131 *
132 */
133 void constructAdjListWithWeights(Graph & g)
134 {
135 // create xadj, adjlist, vwgt, adjwgt and vsize
136 if (Mg.xadj != NULL)
137 {delete [] Mg.xadj;}
138 if (Mg.adjncy != NULL)
139 {delete [] Mg.adjncy;}
140 if (Mg.vwgt != NULL)
141 {delete [] Mg.vwgt;}
142 if (Mg.adjwgt != NULL)
143 {delete [] Mg.adjwgt;}
144 if (Mg.vsize != NULL)
145 {delete [] Mg.vsize;}
146 Mg.xadj = new idx_t[g.getNVertex() + 1];
147 Mg.adjncy = new idx_t[g.getNEdge()];
148 Mg.vwgt = new idx_t[g.getNVertex()];
149 Mg.adjwgt = new idx_t[g.getNEdge()];
150 Mg.vsize = new idx_t[g.getNVertex()];
151
152 //! starting point in the adjacency list
153 size_t prev = 0;
154
155 // actual position
156 size_t id = 0;
157
158 // for each vertex calculate the position of the starting point in the adjacency list
159 for (size_t i = 0; i < g.getNVertex(); i++)
160 {
161 // Add weight to vertex and migration cost
162 Mg.vwgt[i] = g.vertex(i).template get<nm_v_computation>();
163 Mg.vwgt[i] = (Mg.vwgt[i] == 0)?1:Mg.vwgt[i];
164 Mg.vsize[i] = g.vertex(i).template get<nm_v_migration>();
165 Mg.vsize[i] = (Mg.vsize[i] == 0)?1:Mg.vsize[i];
166
167 // Calculate the starting point in the adjacency list
168 Mg.xadj[id] = prev;
169
170 // Create the adjacency list
171 for (size_t s = 0; s < g.getNChilds(i); s++)
172 {
173 Mg.adjncy[prev + s] = g.getChild(i, s);
174
175 // zero values on Metis are dangerous
176 Mg.adjwgt[prev + s] = g.getChildEdge(i, s).template get<nm_e::communication>();
177 Mg.adjwgt[prev + s] = (Mg.adjwgt[prev + s] == 0)?1:Mg.adjwgt[prev + s];
178 }
179
180 // update the position for the next vertex
181 prev += g.getNChilds(i);
182
183 id++;
184 }
185
186 // Fill the last point
187 Mg.xadj[id] = prev;
188 }
189
190
191 void reset_pointer()
192 {
193 Mg.nvtxs = NULL;
194 Mg.ncon = NULL;
195 Mg.xadj = NULL;
196 Mg.adjncy = NULL;
197 Mg.vwgt = NULL;
198 Mg.adjwgt = NULL;
199 Mg.nparts = NULL;
200 Mg.tpwgts = NULL;
201 Mg.ubvec = NULL;
202 Mg.options = NULL;
203 Mg.objval = NULL;
204 Mg.part = NULL;
205 Mg.vsize = NULL;
206 }
207
208public:
209
210 /*! \brief Constructor
211 *
212 * Construct a metis graph from Graph_CSR
213 *
214 * \param g Graph we want to convert to decompose
215 * \param nc number of partitions
216 * \param useWeights tells if weights are used or not
217 *
218 */
219 Metis(Graph & g, size_t nc, bool useWeights)
220 :g(g),n_dec(0)
221 {
222 reset_pointer();
223 initMetisGraph(nc,useWeights);
224 }
225
226 /*! \brief Constructor
227 *
228 * Construct a metis graph from Graph_CSR
229 *
230 * \param g Graph we want to convert to decompose
231 * \param nc number of partitions
232 *
233 */
234 Metis(Graph & g, size_t nc)
235 :g(g),n_dec(0)
236 {
237 reset_pointer();
238 initMetisGraph(nc,false);
239 }
240
241 /*! \brief Constructor
242 *
243 * This constructor does not initialize the internal metis graph
244 * you have to use initMetisGraph to initialize
245 *
246 * \param g Graph we want to convert to decompose
247 *
248 */
249 Metis(Graph & g)
250 :g(g),n_dec(0)
251 {
252 reset_pointer();
253 }
254
255
256 /*! \brief Initialize the METIS graph
257 *
258 * \param nc number of partitions
259 * \param useWeights use the weights on the graph
260 *
261 */
262 void initMetisGraph(int nc, bool useWeights)
263 {
264
265 // Get the number of vertex
266
267 if (Mg.nvtxs != NULL)
268 {delete[] Mg.nvtxs;}
269 Mg.nvtxs = new idx_t[1];
270 Mg.nvtxs[0] = g.getNVertex();
271
272 // Set the number of constrains
273
274 if (Mg.ncon != NULL)
275 {delete[] Mg.ncon;}
276 Mg.ncon = new idx_t[1];
277 Mg.ncon[0] = 1;
278
279 // Set to null the weight of the vertex
280
281 if (Mg.vwgt != NULL)
282 {delete[] Mg.vwgt;}
283 Mg.vwgt = NULL;
284
285 // Put the total communication size to NULL
286
287 if (Mg.vsize != NULL)
288 {delete[] Mg.vsize;}
289 Mg.vsize = NULL;
290
291 // Set to null the weight of the edge
292
293 if (Mg.adjwgt != NULL)
294 {delete[] Mg.adjwgt;}
295 Mg.adjwgt = NULL;
296
297 // construct the adjacency list
298 if (useWeights)
299 constructAdjListWithWeights(g);
300 else
301 constructAdjList(g);
302
303 // Set the total number of partitions
304
305 if (Mg.nparts != NULL)
306 {delete[] Mg.nparts;}
307 Mg.nparts = new idx_t[1];
308 Mg.nparts[0] = nc;
309
310 // Set to null the desired weight for each partition (one for each constrain)
311
312 if (Mg.tpwgts != NULL)
313 {delete[] Mg.tpwgts;}
314 Mg.tpwgts = NULL;
315
316 //! Set to null the partition load imbalance tolerance
317 if (Mg.ubvec != NULL)
318 {delete[] Mg.ubvec;}
319 Mg.ubvec = NULL;
320
321 //! Set tp null additional option for the graph partitioning
322
323 if (Mg.options != NULL)
324 {delete[] Mg.options;}
325 Mg.options = NULL;
326
327 //! set the objective value
328 if (Mg.objval != NULL)
329 {delete[] Mg.objval;}
330 Mg.objval = new idx_t[1];
331
332 //! Is an output vector containing the partition for each vertex
333 if (Mg.part != NULL)
334 {delete[] Mg.part;}
335 Mg.part = new idx_t[g.getNVertex()];
336
337 for (size_t i = 0; i < g.getNVertex(); i++)
338 Mg.part[i] = 0;
339 }
340
341 /*! \brief destructor
342 *
343 * Destructor, It destroy all the memory allocated
344 *
345 */
346 ~Metis()
347 {
348 // Deallocate the Mg structure
349 if (Mg.nvtxs != NULL)
350 {
351 delete[] Mg.nvtxs;
352 }
353
354 if (Mg.ncon != NULL)
355 {
356 delete[] Mg.ncon;
357 }
358
359 if (Mg.xadj != NULL)
360 {
361 delete[] Mg.xadj;
362 }
363
364 if (Mg.adjncy != NULL)
365 {
366 delete[] Mg.adjncy;
367 }
368 if (Mg.vwgt != NULL)
369 {
370 delete[] Mg.vwgt;
371 }
372 if (Mg.adjwgt != NULL)
373 {
374 delete[] Mg.adjwgt;
375 }
376 if (Mg.nparts != NULL)
377 {
378 delete[] Mg.nparts;
379 }
380 if (Mg.tpwgts != NULL)
381 {
382 delete[] Mg.tpwgts;
383 }
384 if (Mg.ubvec != NULL)
385 {
386 delete[] Mg.ubvec;
387 }
388 if (Mg.options != NULL)
389 {
390 delete[] Mg.options;
391 }
392 if (Mg.objval != NULL)
393 {
394 delete[] Mg.objval;
395 }
396 if (Mg.part != NULL)
397 {
398 delete[] Mg.part;
399 }
400 if (Mg.vsize != NULL)
401 {
402 delete[] Mg.vsize;
403 }
404 }
405
406 /*! \brief Decompose the graph
407 *
408 * \tparam i which property store the decomposition
409 *
410 */
411
412 template<unsigned int i>
413 void decompose()
414 {
415 if (Mg.nparts[0] != 1)
416 {
417 // Decompose
418 METIS_PartGraphRecursive(Mg.nvtxs, Mg.ncon, Mg.xadj, Mg.adjncy, Mg.vwgt, Mg.vsize, Mg.adjwgt, Mg.nparts, Mg.tpwgts, Mg.ubvec, Mg.options, Mg.objval, Mg.part);
419
420 // vertex id
421
422 size_t id = 0;
423
424 // For each vertex store the processor that contain the data
425
426 auto it = g.getVertexIterator();
427
428 while (it.isNext())
429 {
430 g.vertex(it).template get<i>() = Mg.part[id];
431
432 ++id;
433 ++it;
434 }
435 }
436 else
437 {
438 // Trivially assign all the domains to the processor 0
439
440 auto it = g.getVertexIterator();
441
442 while (it.isNext())
443 {
444 g.vertex(it).template get<i>() = 0;
445
446 ++it;
447 }
448 }
449
450 n_dec++;
451 }
452
453 /*! \brief Decompose the graph
454 *
455 * \tparam i which property store the decomposition
456 *
457 */
458
459 template<unsigned int i, typename Graph_part>
460 void decompose(Graph_part & gp)
461 {
462 // Decompose
463 METIS_PartGraphRecursive(Mg.nvtxs, Mg.ncon, Mg.xadj, Mg.adjncy, Mg.vwgt, Mg.vsize, Mg.adjwgt, Mg.nparts, Mg.tpwgts, Mg.ubvec, Mg.options, Mg.objval, Mg.part);
464
465 // vertex id
466
467 size_t id = 0;
468
469 // For each vertex store the processor that contain the data
470
471 auto it = gp.getVertexIterator();
472
473 while (it.isNext())
474 {
475 gp.vertex(it).template get<i>() = Mg.part[id];
476
477 ++id;
478 ++it;
479 }
480
481 n_dec++;
482 }
483
484 /*! \brief It set Metis on test
485 *
486 * \param testing set to true to disable the testing
487 *
488 * At the moment disable the seed randomness to keep the result
489 * reproducible
490 *
491 */
492 void onTest(bool testing)
493 {
494 if (testing == false)
495 return;
496
497 if (Mg.options == NULL)
498 {
499 // allocate
500 Mg.options = new idx_t[METIS_NOPTIONS];
501
502 // set default options
503 METIS_SetDefaultOptions(Mg.options);
504 }
505
506 Mg.options[METIS_OPTION_SEED] = 0;
507 }
508
509 /*! \brief Distribution tolerance
510 *
511 * \param tol tolerance
512 *
513 */
514 void setDistTol(real_t tol)
515 {
516 dist_tol = tol;
517 }
518
519 /*! \brief Get the decomposition counter
520 *
521 * \return the decomposition counter
522 *
523 */
524 size_t get_ndec()
525 {
526 return n_dec;
527 }
528
529 /*! \brief Increment the decomposition counter
530 *
531 *
532 */
533 void inc_dec()
534 {
535 n_dec++;
536 }
537};
538
539#endif
540