1 /*! \file MetisSplitter.cpp
2 \brief to call grpah decomposer : METIS
3 \author Xavier Juvigny, ONERA
4 \date Jul. 2nd 2012
5 \date Nov. 30th 2016
6 \date Oct. 15th 2017
7 */
8
9 // This file is part of Dissection
10 //
11 // Dissection is free software: you can redistribute it and/or modify
12 // it under the terms of the GNU General Public License as published by
13 // the Free Software Foundation, either version 3 of the License, or
14 // (at your option) any later version.
15 //
16 // Linking Dissection statically or dynamically with other modules is making
17 // a combined work based on Disssection. Thus, the terms and conditions of
18 // the GNU General Public License cover the whole combination.
19 //
20 // As a special exception, the copyright holders of Dissection give you
21 // permission to combine Dissection program with free software programs or
22 // libraries that are released under the GNU LGPL and with independent modules
23 // that communicate with Dissection solely through the Dissection-fortran
24 // interface. You may copy and distribute such a system following the terms of
25 // the GNU GPL for Dissection and the licenses of the other code concerned,
26 // provided that you include the source code of that other code when and as
27 // the GNU GPL requires distribution of source code and provided that you do
28 // not modify the Dissection-fortran interface.
29 //
30 // Note that people who make modified versions of Dissection are not obligated
31 // to grant this special exception for their modified versions; it is their
32 // choice whether to do so. The GNU General Public License gives permission to
33 // release a modified version without this exception; this exception also makes
34 // it possible to release a modified version which carries forward this
35 // exception. If you modify the Dissection-fortran interface, this exception
36 // does not apply to your modified version of Dissection, and you must remove
37 // this exception when you distribute your modified version.
38 //
39 // This exception is an additional permission under section 7 of the GNU
40 // General Public License, version 3 ("GPLv3")
41 //
42 // Dissection is distributed in the hope that it will be useful,
43 // but WITHOUT ANY WARRANTY; without even the implied warranty of
44 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
45 // GNU General Public License for more details.
46 //
47 // You should have received a copy of the GNU General Public License
48 // along with Dissection. If not, see <http://www.gnu.org/licenses/>.
49 //
50
51 // ============================================================
52 // == Implementation of splitter using Metis library ==
53 // ============================================================
54 #include <stdexcept>
55 #include <cstdio>
56 #include <cstring>
57 #include <algorithm>
58 #ifndef NO_METIS
59 #include "metis.h"
60 #endif
61
62 #include "Compiler/DebugUtils.hpp"
63 #include "Splitters/BitManipulations.hpp"
64 #include "Splitters/MetisSplitter.hpp"
65
compBegOfDomains(unsigned invLevel,unsigned & begDom,unsigned indDom,const int * sizeOfDomains,int * ptOnDomains)66 static unsigned compBegOfDomains(unsigned invLevel,
67 unsigned&begDom,
68 unsigned indDom,
69 const int* sizeOfDomains,
70 int* ptOnDomains)
71 {
72 if (invLevel==1) {
73 ptOnDomains[indDom-1] = begDom;
74 return sizeOfDomains[indDom-1];
75 }
76 else {
77 begDom += compBegOfDomains(invLevel-1, begDom,
78 2*indDom, sizeOfDomains,
79 ptOnDomains);
80 begDom += compBegOfDomains(invLevel-1, begDom,
81 2*indDom+1, sizeOfDomains,
82 ptOnDomains);
83 ptOnDomains[indDom-1] = begDom;
84 return sizeOfDomains[indDom-1];
85 }
86 }
87
88
89 bool
MetisSplitter(unsigned dim,const int * ptRows,const int * indCols,unsigned & nbMaxLevels,unsigned minSize,int * loc2glob,int * glob2loc,int & nbDoms,int * & ptOnDomains,int * & sizeOfDomains,bool checkData,const bool verbose,FILE * fp)90 MetisSplitter(unsigned dim,
91 const int* ptRows, const int* indCols,
92 unsigned& nbMaxLevels, unsigned minSize,
93 int* loc2glob, int* glob2loc, int& nbDoms,
94 int*& ptOnDomains, int*& sizeOfDomains,
95 bool checkData, const bool verbose, FILE *fp)
96 {
97 #ifdef NO_METIS
98 if (verbose) {
99 fprintf(stderr, "%s %d : Metis is not linked\n", __FILE__, __LINE__);
100 }
101 return false;
102 #else
103 /** Check inputs */
104 CHECK(ptRows !=NULL, "Null pointer for ptRows !");
105 CHECK(indCols !=NULL, "Null pointer for indCols !");
106 CHECK(loc2glob!=NULL, "Null pointer for loc2glob!");
107 CHECK(glob2loc!=NULL, "Null pointer for glob2loc!");
108 int ierr;
109 // Compute nbMaxLevels according to the minSize parameter
110 int maxBlocks = dim/minSize;
111 unsigned nbLevels = highestbit(maxBlocks)-1;
112 nbMaxLevels = (nbMaxLevels<nbLevels ? nbMaxLevels : nbLevels);
113 nbDoms = 1<<(nbMaxLevels-1);
114 if (dim==0) return true;
115 idx_t options[METIS_NOPTIONS];
116 METIS_SetDefaultOptions(options);
117 options[METIS_OPTION_PTYPE ] = METIS_PTYPE_RB;
118 options[METIS_OPTION_OBJTYPE ] = METIS_OBJTYPE_NODE;
119 options[METIS_OPTION_IPTYPE ] = METIS_IPTYPE_NODE;
120 options[METIS_OPTION_CTYPE ] = METIS_CTYPE_RM;
121 options[METIS_OPTION_RTYPE ] = METIS_RTYPE_SEP1SIDED;
122 options[METIS_OPTION_NSEPS ] = 1;/*Default value */
123 options[METIS_OPTION_NITER ] = 10;/* Default value */
124 options[METIS_OPTION_UFACTOR ] = 1;
125 options[METIS_OPTION_COMPRESS ] = 0;/* No compress */
126 options[METIS_OPTION_CCORDER ] = 0;/* No connect component detection */
127 options[METIS_OPTION_SEED ] = -1;/* Seed of random algo */
128 options[METIS_OPTION_PFACTOR ] = 0;/* No dense nodes removal */
129 options[METIS_OPTION_NUMBERING] = (idx_t)0; //offset;
130 # ifdef DEBUG
131 options[METIS_OPTION_DBGLVL ] = 255;/* Full debug */
132 # else
133 if(verbose) {
134 options[METIS_OPTION_DBGLVL ] = METIS_DBG_TIME; /* only time profiling */
135 }
136 else{
137 options[METIS_OPTION_DBGLVL ] = 0;
138 }
139 # endif
140 // Opt [0] = 0 ; /* use defaults */
141 // Opt [1] = 3 ; /* matching type */
142 // Opt [2] = 1 ; /* init. partitioning algo*/
143 // Opt [3] = 2 ; /* refinement algorithm */
144 // Opt [4] = 0 ; /* no debug */
145 // Opt [5] = 0 ; /* unused */
146 // Opt [6] = 0 ; /* unused */
147 // Opt [7] = -1 ; /* random seed */
148 //
149 idx_t nvtxs = dim;
150 idx_t *xadj, *adjncy, *perm, *iperm;
151 if (sizeof(idx_t)==sizeof(int)) {// Ok, simple pointer may be ok
152 xadj = (idx_t*)ptRows;
153 adjncy = (idx_t*)indCols;
154 perm = (idx_t*)loc2glob;
155 iperm = (idx_t*)glob2loc;
156 }
157 else {// Arf, some copy to do :(
158 xadj = new idx_t[dim+1];
159 for (unsigned ix = 0; ix <= dim; ix++)
160 xadj[ix] = ptRows[ix];
161 adjncy = new idx_t[ptRows[dim]];
162 for (unsigned ix = 0; ix < ptRows[dim]; ++ix)
163 adjncy[ix] = indCols[ix];
164 perm = new idx_t[dim];
165 iperm = new idx_t[dim];
166 }
167 if (verbose) {
168 fprintf(fp, "%s %d : nbDoms = %d\n",
169 __FILE__, __LINE__, nbDoms);
170 }
171 idx_t* sizes = new idx_t[2*nbDoms-1];
172 /** This metis subroutine make dissection with some history
173 in sizes array
174 */
175 int err = METIS_NodeNDP(nvtxs, xadj, adjncy, NULL,
176 idx_t(nbDoms), options,
177 iperm, perm, sizes);
178 if (verbose && (err != METIS_OK)) {
179 fprintf(fp, "Failed to complete Metis bisection !\n");
180 if (err == METIS_ERROR_INPUT)
181 fprintf(fp, "\t Wrong arguments given to Metis\n");
182 else if (err == METIS_ERROR_MEMORY)
183 fprintf(fp, "\t Metis fails allocate some memory\n");
184 else
185 fprintf(fp, "\t Other error in Metis\n");
186 }
187 // TRACE("Finish to split ;)\n");
188
189 nbDoms = nbDoms*2-1;
190 // Compute where start each domain per bisection level
191 // and compute the size of each subdomain
192 ptOnDomains = new int[nbDoms+1];
193 sizeOfDomains = new int[nbDoms];
194 unsigned ptSize = 0;
195 unsigned begDom = 0;
196 // Metis store size of each domains and interface
197 // per level of tree, beginning with the leaves of the tree
198 for (unsigned iLvl = nbMaxLevels; iLvl > 0; iLvl --)
199 {
200 for (unsigned iDom = (1<<(iLvl-1)); iDom < (1<<iLvl);
201 ++iDom)
202 {
203 assert(ptSize < nbDoms);
204 sizeOfDomains[iDom-1] = sizes[ptSize];
205 ptSize++;
206 }
207 }
208 delete [] sizes;
209 // But renumbering is done as left child - right child - father
210 // way !
211 compBegOfDomains(nbMaxLevels, begDom, 1, sizeOfDomains,ptOnDomains);
212 ptOnDomains[nbDoms] = dim;
213 //assert(begDom==dim);
214 assert(ptOnDomains[0]+sizeOfDomains[0]==dim);
215 if (sizeof(idx_t)!=sizeof(int)) {
216 for (unsigned i = 0; i < dim; i++)
217 {
218 loc2glob[i] = perm[i];
219 glob2loc[i] = iperm[i];
220 }
221 delete [] iperm;
222 delete [] perm;
223 delete [] adjncy;
224 delete [] xadj;
225 }
226 return true;
227 #endif
228 }
229