45 #ifndef _ZOLTAN2_ALGSCOTCH_HPP_ 46 #define _ZOLTAN2_ALGSCOTCH_HPP_ 59 #ifndef HAVE_ZOLTAN2_SCOTCH 65 template <
typename Adapter>
70 const RCP<
const Comm<int> > &problemComm,
74 throw std::runtime_error(
75 "BUILD ERROR: Scotch requested but not compiled into Zoltan2.\n" 76 "Please set CMake flag Zoltan2_ENABLE_Scotch:BOOL=ON.");
85 #ifdef HAVE_ZOLTAN2_SCOTCH 92 #ifndef HAVE_ZOLTAN2_MPI 99 #ifdef SHOW_ZOLTAN2_SCOTCH_MEMORY 112 #ifdef HAVE_SCOTCH_GETMEMORYMAX 114 extern size_t SCOTCH_getMemoryMax();
117 #error "Either turn off ZOLTAN2_ENABLE_SCOTCH_MEMORY_REPORT in cmake configure, or see SHOW_ZOLTAN2_SCOTCH_MEMORY comment in Zoltan2_AlgScotch.hpp" 118 #endif // HAVE_SCOTCH_GETMEMORYMAX 119 #endif // SHOW_ZOLTAN2_SCOTCH_MEMORY 121 template <
typename Adapter>
127 typedef typename Adapter::lno_t
lno_t;
128 typedef typename Adapter::gno_t
gno_t;
129 typedef typename Adapter::scalar_t
scalar_t;
130 typedef typename Adapter::part_t
part_t;
143 const RCP<
const Comm<int> > &problemComm__,
144 const RCP<graphModel_t> &model__) :
145 env(env__), problemComm(problemComm__),
146 #ifdef HAVE_ZOLTAN2_MPI 147 mpicomm(Teuchos::getRawMpiComm(*problemComm__)),
156 const RCP<const Environment> env;
157 const RCP<const Comm<int> > problemComm;
158 #ifdef HAVE_ZOLTAN2_MPI 159 const MPI_Comm mpicomm;
161 const RCP<GraphModel<typename Adapter::base_adapter_t> > model;
169 template <
typename Adapter>
176 size_t numGlobalParts = solution->getTargetGlobalNumberOfParts();
178 SCOTCH_Num partnbr=0;
181 #ifdef HAVE_ZOLTAN2_MPI 183 int me = problemComm->getRank();
185 const SCOTCH_Num baseval = 0;
188 SCOTCH_Strat stratstr;
190 SCOTCH_stratInit(&stratstr);
193 SCOTCH_Dgraph *gr = SCOTCH_dgraphAlloc();
194 ierr = SCOTCH_dgraphInit(gr, mpicomm);
196 env->globalInputAssertion(__FILE__, __LINE__,
"SCOTCH_dgraphInit",
200 ArrayView<const gno_t> vtxID;
201 ArrayView<StridedData<lno_t, scalar_t> > vwgts;
202 size_t nVtx = model->getVertexList(vtxID, vwgts);
203 SCOTCH_Num vertlocnbr=0;
205 SCOTCH_Num vertlocmax = vertlocnbr;
208 ArrayView<const gno_t> edgeIds;
209 ArrayView<const lno_t> offsets;
210 ArrayView<StridedData<lno_t, scalar_t> > ewgts;
212 size_t nEdge = model->getEdgeList(edgeIds, offsets, ewgts);
214 SCOTCH_Num edgelocnbr=0;
216 const SCOTCH_Num edgelocsize = edgelocnbr;
218 SCOTCH_Num *vertloctab;
221 SCOTCH_Num *edgeloctab;
225 SCOTCH_Num *vendloctab = NULL;
226 SCOTCH_Num *vlblloctab = NULL;
227 SCOTCH_Num *edgegsttab = NULL;
230 SCOTCH_Num *velotab = NULL;
231 SCOTCH_Num *edlotab = NULL;
233 int nVwgts = model->getNumWeightsPerVertex();
234 int nEwgts = model->getNumWeightsPerEdge();
235 if (nVwgts > 1 && me == 0) {
236 std::cerr <<
"Warning: NumWeightsPerVertex is " << nVwgts
237 <<
" but Scotch allows only one weight. " 238 <<
" Zoltan2 will use only the first weight per vertex." 241 if (nEwgts > 1 && me == 0) {
242 std::cerr <<
"Warning: NumWeightsPerEdge is " << nEwgts
243 <<
" but Scotch allows only one weight. " 244 <<
" Zoltan2 will use only the first weight per edge." 249 velotab =
new SCOTCH_Num[nVtx+1];
251 scale_weights(nVtx, vwgts[0], velotab);
255 edlotab =
new SCOTCH_Num[nEdge+1];
257 scale_weights(nEdge, ewgts[0], edlotab);
261 ierr = SCOTCH_dgraphBuild(gr, baseval, vertlocnbr, vertlocmax,
262 vertloctab, vendloctab, velotab, vlblloctab,
263 edgelocnbr, edgelocsize,
264 edgeloctab, edgegsttab, edlotab);
266 env->globalInputAssertion(__FILE__, __LINE__,
"SCOTCH_dgraphBuild",
270 ArrayRCP<part_t> partList(
new part_t[nVtx], 0, nVtx,
true);
271 SCOTCH_Num *partloctab = NULL;
272 if (nVtx && (
sizeof(SCOTCH_Num) ==
sizeof(
part_t))) {
274 partloctab = (SCOTCH_Num *) partList.getRawPtr();
280 partloctab =
new SCOTCH_Num[nVtx+1];
284 float *partsizes =
new float[numGlobalParts];
285 if (!solution->criteriaHasUniformPartSizes(0))
286 for (
size_t i=0; i<numGlobalParts; i++)
287 partsizes[i] = solution->getCriteriaPartSize(0, i);
289 for (
size_t i=0; i<numGlobalParts; i++)
290 partsizes[i] = 1.0 /
float(numGlobalParts);
294 SCOTCH_archInit(&archdat);
296 SCOTCH_Num velosum = 0;
297 SCOTCH_dgraphSize (gr, &velosum, NULL, NULL, NULL);
298 SCOTCH_Num *goalsizes =
new SCOTCH_Num[partnbr];
304 for (SCOTCH_Num i = 0; i < partnbr; i++)
305 goalsizes[i] = SCOTCH_Num(ceil(velosum * partsizes[i]));
308 SCOTCH_archCmpltw(&archdat, partnbr, goalsizes);
311 ierr = SCOTCH_dgraphMap(gr, &archdat, &stratstr, partloctab);
313 env->globalInputAssertion(__FILE__, __LINE__,
"SCOTCH_dgraphMap",
316 SCOTCH_archExit(&archdat);
321 #ifdef SHOW_ZOLTAN2_SCOTCH_MEMORY 322 int me = env->comm_->getRank();
325 #ifdef HAVE_SCOTCH_ZOLTAN2_GETMEMORYMAX 327 size_t scotchBytes = SCOTCH_getMemoryMax();
328 std::cout <<
"Rank " << me <<
": Maximum bytes used by Scotch: ";
329 std::cout << scotchBytes << std::endl;
334 SCOTCH_dgraphExit(gr);
336 SCOTCH_stratExit(&stratstr);
340 if ((
sizeof(SCOTCH_Num) !=
sizeof(
part_t)) || (nVtx == 0)) {
341 for (
size_t i = 0; i < nVtx; i++) partList[i] = partloctab[i];
342 delete [] partloctab;
345 solution->setParts(partList);
347 env->memory(
"Zoltan2-Scotch: After creating solution");
353 if (nVwgts)
delete [] velotab;
354 if (nEwgts)
delete [] edlotab;
356 #else // DO NOT HAVE MPI 364 ArrayView<const gno_t> vtxID;
365 ArrayView<StridedData<lno_t, scalar_t> > vwgts;
366 size_t nVtx = model->getVertexList(vtxID, vwgts);
368 ArrayRCP<part_t> partList(
new part_t[nVtx], 0, nVtx,
true);
369 for (
size_t i = 0; i < nVtx; i++) partList[i] = 0;
371 solution->setParts(partList);
373 #endif // DO NOT HAVE MPI 386 template <
typename Adapter>
393 const double INT_EPSILON = 1e-5;
395 SCOTCH_Num nonint, nonint_local = 0;
396 double sum_wgt, sum_wgt_local = 0.;
397 double max_wgt, max_wgt_local = 0.;
401 for (
size_t i = 0; i < n; i++) {
402 double fw = double(fwgts[i]);
404 SCOTCH_Num tmp = (SCOTCH_Num) floor(fw + .5);
405 if (fabs((
double)tmp-fw) > INT_EPSILON) {
410 if (fw > max_wgt_local) max_wgt_local = fw;
413 Teuchos::reduceAll<int,SCOTCH_Num>(*problemComm, Teuchos::REDUCE_MAX, 1,
414 &nonint_local, &nonint);
415 Teuchos::reduceAll<int,double>(*problemComm, Teuchos::REDUCE_SUM, 1,
416 &sum_wgt_local, &sum_wgt);
417 Teuchos::reduceAll<int,double>(*problemComm, Teuchos::REDUCE_MAX, 1,
418 &max_wgt_local, &max_wgt);
421 const double max_wgt_sum = double(SCOTCH_NUMMAX/8);
425 if (nonint || (max_wgt <= INT_EPSILON) || (sum_wgt > max_wgt_sum)) {
427 if (sum_wgt != 0.) scale = max_wgt_sum/sum_wgt;
431 for (
size_t i = 0; i < n; i++)
432 iwgts[i] = (SCOTCH_Num) ceil(
double(fwgts[i])*scale);
438 #endif // HAVE_ZOLTAN2_SCOTCH
static void DELETE_TPL_T_ARRAY(tpl_t **a)
fast typical checks for valid arguments
Defines the PartitioningSolution class.
AlgPTScotch(const RCP< const Environment > &env, const RCP< const Comm< int > > &problemComm, const RCP< GraphModel< typename Adapter::base_adapter_t > > &model)
Adapter::scalar_t scalar_t
A PartitioningSolution is a solution to a partitioning problem.
The StridedData class manages lists of weights or coordinates.
Algorithm defines the base class for all algorithms.
static void ASSIGN_TPL_T_ARRAY(tpl_t **a, ArrayView< zno_t > &b)
virtual void partition(const RCP< PartitioningSolution< Adapter > > &solution)
Partitioning method.
Traits class to handle conversions between gno_t/lno_t and TPL data types (e.g., ParMETIS's idx_t...
GraphModel defines the interface required for graph models.
static void ASSIGN_TPL_T(tpl_t &a, zno_t b)
Defines the GraphModel interface.
A gathering of useful namespace methods.