51 #ifndef ZOLTAN2_METRIC_HPP 52 #define ZOLTAN2_METRIC_HPP 56 #include <Tpetra_Map.hpp> 57 #include <Tpetra_CrsMatrix.hpp> 74 template <
typename scalar_t>
81 values_ = arcp(tmp, 0, evalNumMetrics,
true);
83 ArrayRCP<scalar_t> values_;
84 std::string metricName_;
131 void setName(std::string name) { metricName_ = name;}
164 const std::string &
getName()
const {
return metricName_; }
200 template <
typename scalar_t>
207 values_ = arcp(tmp, 0, evalNumMetrics,
true);
209 ArrayRCP<scalar_t> values_;
210 std::string metricName_;
231 values_(), metricName_(mname) {
236 values_(), metricName_(
"unset") {
240 void setName(std::string name) { metricName_ = name;}
249 const std::string &
getName()
const {
return metricName_; }
259 template <
typename scalar_t>
262 std::string label(metricName_);
265 std::ostringstream oss;
268 oss << metricName_ <<
" (1)";
271 oss << metricName_ <<
" (2)";
274 oss << metricName_ <<
" (inf)";
277 oss << metricName_ <<
" (?)";
284 os << std::setw(20) << label;
285 os << std::setw(12) << std::setprecision(4) << values_[
evalGlobalMin];
286 os << std::setw(12) << std::setprecision(4) << values_[
evalGlobalMax];
287 os << std::setw(12) << std::setprecision(4) << values_[
evalGlobalAvg];
288 os << std::setw(2) <<
" ";
295 template <
typename scalar_t>
298 os << std::setw(20) <<
" ";
299 os << std::setw(36) <<
"------------SUM PER PART-----------";
300 os << std::setw(2) <<
" ";
301 os << std::setw(18) <<
"IMBALANCE PER PART";
304 os << std::setw(20) <<
" ";
305 os << std::setw(12) <<
"min" << std::setw(12) <<
"max" << std::setw(12) <<
"avg";
306 os << std::setw(2) <<
" ";
307 os << std::setw(6) <<
"min" << std::setw(6) <<
"max" << std::setw(6) <<
"avg";
311 template <
typename scalar_t>
314 std::string label(metricName_);
316 os << std::setw(20) << label;
317 os << std::setw(12) << std::setprecision(4) << values_[
evalGlobalSum];
318 os << std::setw(12) << std::setprecision(4) << values_[
evalGlobalMax];
322 template <
typename scalar_t>
325 os << std::setw(20) <<
" ";
326 os << std::setw(24) <<
"----------SUM----------";
329 os << std::setw(20) <<
" ";
330 os << std::setw(12) <<
"sum" << std::setw(12) <<
"max";
347 template <
typename scalar_t>
349 int offset, scalar_t &min, scalar_t &max, scalar_t &sum)
351 if (v.size() < 1)
return;
352 min = max = sum = v[offset];
354 for (
int i=offset+stride; i < v.size(); i += stride){
355 if (v[i] < min) min = v[i];
356 else if (v[i] > max) max = v[i];
369 template <
typename scalar_t>
371 int offset, scalar_t &max, scalar_t &sum)
373 if (v.size() < 1)
return;
374 max = sum = v[offset];
376 for (
int i=offset+stride; i < v.size(); i += stride){
377 if (v[i] > max) max = v[i];
400 template <
typename scalar_t,
typename lno_t,
typename part_t>
402 const RCP<const Environment> &env,
403 part_t numberOfParts,
404 const ArrayView<const part_t> &parts,
409 env->localInputAssertion(__FILE__, __LINE__,
"parts or weights",
412 int numObjects = parts.size();
413 int vwgtDim = vwgts.size();
415 memset(weights, 0,
sizeof(scalar_t) * numberOfParts);
421 for (lno_t i=0; i < numObjects; i++){
425 else if (vwgtDim == 1){
426 for (lno_t i=0; i < numObjects; i++){
427 weights[parts[i]] += vwgts[0][i];
433 for (
int wdim=0; wdim < vwgtDim; wdim++){
434 for (lno_t i=0; i < numObjects; i++){
435 weights[parts[i]] += vwgts[wdim][i];
441 for (lno_t i=0; i < numObjects; i++){
443 for (
int wdim=0; wdim < vwgtDim; wdim++)
444 ssum += (vwgts[wdim][i] * vwgts[wdim][i]);
445 weights[parts[i]] += sqrt(ssum);
450 for (lno_t i=0; i < numObjects; i++){
452 for (
int wdim=0; wdim < vwgtDim; wdim++)
453 if (vwgts[wdim][i] > max)
454 max = vwgts[wdim][i];
455 weights[parts[i]] += max;
460 env->localBugAssertion(__FILE__, __LINE__,
"invalid norm",
false,
508 template <
typename scalar_t,
typename lno_t,
typename part_t>
510 const RCP<const Environment> &env,
511 const RCP<
const Comm<int> > &comm,
512 const ArrayView<const part_t> &part,
517 part_t &numNonemptyParts,
519 ArrayRCP<scalar_t> &globalSums)
525 numParts = numNonemptyParts = 0;
528 if (vwgtDim) numMetrics++;
529 if (vwgtDim > 1) numMetrics += vwgtDim;
532 mv_t *newMetrics =
new mv_t [numMetrics];
533 env->localMemoryAssertion(__FILE__, __LINE__, numMetrics, newMetrics);
534 ArrayRCP<mv_t> metricArray(newMetrics, 0, numMetrics,
true);
536 metrics = metricArray;
542 lno_t localNumObj = part.size();
543 part_t localNum[2], globalNum[2];
544 localNum[0] =
static_cast<part_t
>(vwgtDim);
547 for (lno_t i=0; i < localNumObj; i++)
548 if (part[i] > localNum[1]) localNum[1] = part[i];
551 reduceAll<int, part_t>(*comm, Teuchos::REDUCE_MAX, 2,
552 localNum, globalNum);
556 env->globalBugAssertion(__FILE__,__LINE__,
557 "inconsistent number of vertex weights",
560 part_t nparts = globalNum[1] + 1;
562 part_t globalSumSize = nparts * numMetrics;
563 scalar_t * sumBuf =
new scalar_t [globalSumSize];
564 env->localMemoryAssertion(__FILE__, __LINE__, globalSumSize, sumBuf);
565 globalSums = arcp(sumBuf, 0, globalSumSize);
570 scalar_t *localBuf =
new scalar_t [globalSumSize];
571 env->localMemoryAssertion(__FILE__, __LINE__, globalSumSize, localBuf);
572 memset(localBuf, 0,
sizeof(scalar_t) * globalSumSize);
574 scalar_t *obj = localBuf;
576 for (lno_t i=0; i < localNumObj; i++)
581 scalar_t *wgt = localBuf + nparts;
583 normedPartWeights<scalar_t, lno_t, part_t>(env, nparts,
584 part, vwgts, mcNorm, wgt);
592 for (
int vdim = 0; vdim < vwgtDim; vdim++){
593 for (lno_t i=0; i < localNumObj; i++)
594 wgt[part[i]] += vwgts[vdim][i];
604 metrics[next].setName(
"object count");
605 metrics[next].setLocalSum(localNumObj);
609 scalar_t *wgt = localBuf + nparts;
610 scalar_t total = 0.0;
612 for (
int p=0; p < nparts; p++){
617 metrics[next].setName(
"weight 1");
619 metrics[next].setName(
"normed weight");
621 metrics[next].setLocalSum(total);
625 for (
int vdim = 0; vdim < vwgtDim; vdim++){
628 for (
int p=0; p < nparts; p++){
632 std::ostringstream oss;
633 oss <<
"weight " << vdim+1;
635 metrics[next].setName(oss.str());
636 metrics[next].setLocalSum(total);
646 reduceAll<int, scalar_t>(*comm, Teuchos::REDUCE_SUM, globalSumSize,
657 scalar_t min=0, max=0, sum=0;
660 ArrayView<scalar_t> objVec(obj, nparts);
661 getStridedStats<scalar_t>(objVec, 1, 0, min, max, sum);
663 metrics[next].setGlobalMin(min);
664 metrics[next].setGlobalMax(max);
665 metrics[next].setGlobalSum(sum);
666 metrics[next].setGlobalAvg(sum / nparts);
670 scalar_t *wgt = sumBuf + nparts;
672 ArrayView<scalar_t> normedWVec(wgt, nparts);
673 getStridedStats<scalar_t>(normedWVec, 1, 0, min, max, sum);
678 metrics[next].setGlobalMin(min);
679 metrics[next].setGlobalMax(max);
680 metrics[next].setGlobalSum(sum);
681 metrics[next].setGlobalAvg(sum / nparts);
685 for (
int vdim=0; vdim < vwgtDim; vdim++){
687 ArrayView<scalar_t> fromVec(wgt, nparts);
688 getStridedStats<scalar_t>(fromVec, 1, 0, min, max, sum);
690 metrics[next].setGlobalMin(min);
691 metrics[next].setGlobalMax(max);
692 metrics[next].setGlobalSum(sum);
693 metrics[next].setGlobalAvg(sum / nparts);
710 numNonemptyParts = numParts;
712 for (part_t p=0; p < numParts; p++)
713 if (obj[p] == 0) numNonemptyParts--;
747 template <
typename Adapter>
749 const RCP<const Environment> &env,
750 const RCP<
const Comm<int> > &comm,
752 const ArrayView<const typename Adapter::part_t> &part,
753 typename Adapter::part_t &numParts,
755 ArrayRCP<typename Adapter::scalar_t> &globalSums)
763 int ewgtDim = graph->getNumWeightsPerEdge();
766 if (ewgtDim > 1) numMetrics = ewgtDim;
768 typedef typename Adapter::scalar_t scalar_t;
769 typedef typename Adapter::gno_t gno_t;
770 typedef typename Adapter::lno_t lno_t;
771 typedef typename Adapter::node_t node_t;
772 typedef typename Adapter::part_t part_t;
776 typedef Tpetra::CrsMatrix<part_t,lno_t,gno_t,node_t> sparse_matrix_type;
777 typedef Tpetra::Vector<part_t,lno_t,gno_t,node_t> vector_t;
778 typedef Tpetra::Map<lno_t, gno_t, node_t> map_type;
780 const GST
INVALID = Teuchos::OrdinalTraits<GST>::invalid ();
784 mv_t *newMetrics =
new mv_t [numMetrics];
785 env->localMemoryAssertion(__FILE__,__LINE__,numMetrics,newMetrics);
786 ArrayRCP<mv_t> metricArray(newMetrics, 0, numMetrics,
true);
788 metrics = metricArray;
794 lno_t localNumObj = part.size();
795 part_t localNum[2], globalNum[2];
796 localNum[0] =
static_cast<part_t
>(ewgtDim);
799 for (lno_t i=0; i < localNumObj; i++)
800 if (part[i] > localNum[1]) localNum[1] = part[i];
803 reduceAll<int, part_t>(*comm, Teuchos::REDUCE_MAX, 2,
804 localNum, globalNum);
808 env->globalBugAssertion(__FILE__,__LINE__,
809 "inconsistent number of edge weights",
812 part_t nparts = globalNum[1] + 1;
814 part_t globalSumSize = nparts * numMetrics;
815 scalar_t * sumBuf =
new scalar_t [globalSumSize];
816 env->localMemoryAssertion(__FILE__, __LINE__, globalSumSize, sumBuf);
817 globalSums = arcp(sumBuf, 0, globalSumSize);
822 scalar_t *localBuf =
new scalar_t [globalSumSize];
823 env->localMemoryAssertion(__FILE__,__LINE__,globalSumSize,localBuf);
824 memset(localBuf, 0,
sizeof(scalar_t) * globalSumSize);
826 scalar_t *cut = localBuf;
828 ArrayView<const gno_t> Ids;
829 ArrayView<input_t> vwgts;
831 graph->getVertexList(Ids, vwgts);
833 ArrayView<const gno_t> edgeIds;
834 ArrayView<const lno_t> offsets;
835 ArrayView<input_t> wgts;
837 graph->getEdgeList(edgeIds, offsets, wgts);
842 RCP<const map_type> vertexMapG;
845 gno_t min = std::numeric_limits<gno_t>::max();
847 for (lno_t i = 0; i < localNumObj; ++i) {
848 if (Ids[i] < min) min = Ids[i];
849 size_t ncols = offsets[i+1] - offsets[i];
850 if (ncols > maxcols) maxcols = ncols;
854 Teuchos::reduceAll<int, gno_t>(*comm,Teuchos::REDUCE_MIN,1,&min,&gmin);
857 vertexMapG = rcp(
new map_type(INVALID, Ids, gmin, comm));
863 RCP<sparse_matrix_type> adjsMatrix;
866 adjsMatrix = rcp (
new sparse_matrix_type (vertexMapG, 0));
868 Array<part_t> justOneA(maxcols, 1);
870 for (lno_t localElement=0; localElement<localNumObj; ++localElement){
872 size_t ncols = offsets[localElement+1] - offsets[localElement];
873 adjsMatrix->insertGlobalValues(Ids[localElement],
874 edgeIds(offsets[localElement], ncols),
879 adjsMatrix->fillComplete ();
882 RCP<vector_t> scaleVec = Teuchos::rcp(
new vector_t(vertexMapG,
false) );
883 for (lno_t localElement=0; localElement<localNumObj; ++localElement) {
884 scaleVec->replaceLocalValue(localElement,part[localElement]);
888 adjsMatrix->rightScale(*scaleVec);
889 Array<gno_t> Indices;
890 Array<part_t> Values;
893 for (lno_t i=0; i < localNumObj; i++) {
894 const gno_t globalRow = Ids[i];
895 size_t NumEntries = adjsMatrix->getNumEntriesInGlobalRow (globalRow);
896 Indices.resize (NumEntries);
897 Values.resize (NumEntries);
898 adjsMatrix->getGlobalRowCopy (globalRow,Indices(),Values(),NumEntries);
900 for (
size_t j=0; j < NumEntries; j++)
901 if (part[i] != Values[j])
908 scalar_t *wgt = localBuf;
909 for (
int edim = 0; edim < ewgtDim; edim++){
910 for (lno_t i=0; i < localNumObj; i++) {
911 const gno_t globalRow = Ids[i];
912 size_t NumEntries = adjsMatrix->getNumEntriesInGlobalRow (globalRow);
913 Indices.resize (NumEntries);
914 Values.resize (NumEntries);
915 adjsMatrix->getGlobalRowCopy (globalRow,Indices(),Values(),NumEntries);
917 for (
size_t j=0; j < NumEntries; j++)
918 if (part[i] != Values[j])
919 wgt[part[i]] += wgts[edim][offsets[i] + j];
929 reduceAll<int, scalar_t>(*comm, Teuchos::REDUCE_SUM, globalSumSize,
940 scalar_t max=0, sum=0;
943 ArrayView<scalar_t> cutVec(cut, nparts);
944 getStridedStats<scalar_t>(cutVec, 1, 0, max, sum);
946 metrics[next].setName(
"cut count");
947 metrics[next].setGlobalMax(max);
948 metrics[next].setGlobalSum(sum);
951 scalar_t *wgt = sumBuf;
953 for (
int edim=0; edim < ewgtDim; edim++){
954 ArrayView<scalar_t> fromVec(wgt, nparts);
955 getStridedStats<scalar_t>(fromVec, 1, 0, max, sum);
957 std::ostringstream oss;
958 oss <<
"weight " << edim+1;
960 metrics[next].setName(oss.str());
961 metrics[next].setGlobalMax(max);
962 metrics[next].setGlobalSum(sum);
1003 template <
typename scalar_t,
typename part_t>
1005 const scalar_t *psizes, scalar_t sumVals ,
const scalar_t *
vals,
1006 scalar_t &min, scalar_t &max, scalar_t &avg)
1011 if (sumVals <= 0 || targetNumParts < 1 || numParts < 1)
1014 if (targetNumParts==1 || numParts==1){
1015 min = max = avg = 0;
1020 scalar_t target = sumVals / targetNumParts;
1021 for (part_t p=0; p < numParts; p++){
1022 scalar_t diff = vals[p] - target;
1023 scalar_t adiff = (diff >= 0 ? diff : -diff);
1024 scalar_t tmp = adiff / target;
1026 if (tmp > max) max = tmp;
1027 if (tmp < min) min = tmp;
1029 part_t emptyParts = targetNumParts - numParts;
1030 if (emptyParts > 0){
1037 for (part_t p=0; p < targetNumParts; p++){
1040 scalar_t target = sumVals * psizes[p];
1041 scalar_t diff = vals[p] - target;
1042 scalar_t adiff = (diff >= 0 ? diff : -diff);
1043 scalar_t tmp = adiff / target;
1045 if (tmp > max) max = tmp;
1046 if (tmp < min) min = tmp;
1057 avg /= targetNumParts;
1093 template <
typename scalar_t,
typename part_t>
1095 int numSizes, ArrayView<ArrayRCP<scalar_t> > psizes,
1096 scalar_t sumVals ,
const scalar_t *
vals,
1097 scalar_t &min, scalar_t &max, scalar_t &avg)
1102 if (sumVals <= 0 || targetNumParts < 1 || numParts < 1)
1105 if (targetNumParts==1 && numParts==1){
1106 min = max = avg = 0;
1110 bool allUniformParts =
true;
1111 for (
int i=0; i < numSizes; i++){
1112 if (psizes[i].size() != 0){
1113 allUniformParts =
false;
1118 if (allUniformParts){
1119 computeImbalances<scalar_t, part_t>(numParts, targetNumParts, NULL,
1120 sumVals,
vals, min, max, avg);
1124 double uniformSize = 1.0 / targetNumParts;
1125 std::vector<double> sizeVec(numSizes);
1126 for (
int i=0; i < numSizes; i++){
1127 sizeVec[i] = uniformSize;
1130 for (part_t p=0; p < numParts; p++){
1137 if (p >= targetNumParts)
1142 double targetNorm = 0;
1143 for (
int i=0; i < numSizes; i++) {
1144 if (psizes[i].size() > 0)
1145 sizeVec[i] = psizes[i][p];
1146 sizeVec[i] *= sumVals;
1147 targetNorm += (sizeVec[i] * sizeVec[i]);
1149 targetNorm = sqrt(targetNorm);
1154 if (targetNorm > 0){
1158 std::vector<double> actual(numSizes);
1159 double actualNorm = 0.;
1160 for (
int i=0; i < numSizes; i++) {
1161 actual[i] = vals[p] * -1.0;
1162 actual[i] += sizeVec[i];
1163 actualNorm += (actual[i] * actual[i]);
1165 actualNorm = sqrt(actualNorm);
1169 scalar_t imbalance = actualNorm / targetNorm;
1171 if (imbalance < min)
1173 else if (imbalance > max)
1179 part_t numEmptyParts = 0;
1181 for (part_t p=numParts; p < targetNumParts; p++){
1182 bool nonEmptyPart =
false;
1183 for (
int i=0; !nonEmptyPart && i < numSizes; i++)
1184 if (psizes[i].size() > 0 && psizes[i][p] > 0.0)
1185 nonEmptyPart =
true;
1194 if (numEmptyParts > 0){
1195 avg += numEmptyParts;
1200 avg /= targetNumParts;
1235 template <
typename Adapter>
1237 const RCP<const Environment> &env,
1238 const RCP<
const Comm<int> > &comm,
1240 const RCP<const typename Adapter::base_adapter_t> &ia,
1242 bool useDegreeAsWeight,
1244 typename Adapter::part_t &numParts,
1245 typename Adapter::part_t &numNonemptyParts,
1250 typedef typename Adapter::scalar_t scalar_t;
1251 typedef typename Adapter::gno_t gno_t;
1252 typedef typename Adapter::lno_t lno_t;
1253 typedef typename Adapter::part_t part_t;
1254 typedef typename Adapter::base_adapter_t base_adapter_t;
1259 size_t numLocalObjects = ia->getLocalNumIDs();
1263 const part_t *parts;
1264 if (solution != Teuchos::null) {
1265 parts = solution->getPartListView();
1266 env->localInputAssertion(__FILE__, __LINE__,
"parts not set",
1269 part_t *procs =
new part_t [numLocalObjects];
1270 for (
size_t i = 0; i < numLocalObjects; i++) procs[i] = comm->getRank();
1273 ArrayView<const part_t> partArray(parts, numLocalObjects);
1277 int nWeights = ia->getNumWeightsPerID();
1278 int numCriteria = (nWeights > 0 ? nWeights : 1);
1279 Array<sdata_t>
weights(numCriteria);
1284 weights[0] = sdata_t();
1287 if (useDegreeAsWeight) {
1288 ArrayView<const gno_t> Ids;
1289 ArrayView<sdata_t> vwgts;
1290 if (graphModel == Teuchos::null) {
1291 std::bitset<NUM_MODEL_FLAGS> modelFlags;
1292 RCP<GraphModel<base_adapter_t> > graph;
1294 graph->getVertexList(Ids, vwgts);
1296 graphModel->getVertexList(Ids, vwgts);
1298 scalar_t *wgt =
new scalar_t[numLocalObjects];
1299 for (
int i=0; i < nWeights; i++){
1300 for (
size_t j=0; j < numLocalObjects; j++) {
1301 wgt[j] = vwgts[i][j];
1303 ArrayRCP<const scalar_t> wgtArray(wgt,0,numLocalObjects,
false);
1304 weights[i] = sdata_t(wgtArray, 1);
1307 for (
int i=0; i < nWeights; i++){
1308 const scalar_t *wgt;
1310 ia->getWeightsView(wgt, stride, i);
1311 ArrayRCP<const scalar_t> wgtArray(wgt,0,stride*numLocalObjects,
false);
1312 weights[i] = sdata_t(wgtArray, stride);
1319 part_t targetNumParts = solution->getTargetGlobalNumberOfParts();
1320 scalar_t *psizes = NULL;
1322 ArrayRCP<ArrayRCP<scalar_t> > partSizes(numCriteria);
1323 for (
int dim=0; dim < numCriteria; dim++){
1324 if (solution->criteriaHasUniformPartSizes(dim) !=
true){
1325 psizes =
new scalar_t [targetNumParts];
1326 env->localMemoryAssertion(__FILE__, __LINE__, numParts, psizes);
1327 for (part_t i=0; i < targetNumParts; i++){
1328 psizes[i] = solution->getCriteriaPartSize(dim, i);
1330 partSizes[dim] = arcp(psizes, 0, targetNumParts,
true);
1338 ArrayRCP<scalar_t> globalSums;
1341 globalSumsByPart<scalar_t, lno_t, part_t>(env, comm,
1342 partArray, nWeights, weights.view(0, numCriteria), mcNorm,
1343 numParts, numNonemptyParts, metrics, globalSums);
1351 scalar_t *objCount = globalSums.getRawPtr();
1352 scalar_t min, max, avg;
1355 if (partSizes[0].size() > 0)
1356 psizes = partSizes[0].getRawPtr();
1358 computeImbalances<scalar_t, part_t>(numParts, targetNumParts, psizes,
1359 metrics[0].getGlobalSum(), objCount,
1362 metrics[0].setMinImbalance(1.0 + min);
1363 metrics[0].setMaxImbalance(1.0 + max);
1364 metrics[0].setAvgImbalance(1.0 + avg);
1369 scalar_t *wgts = globalSums.getRawPtr() + numParts;
1371 if (metrics.size() > 1){
1373 computeImbalances<scalar_t, part_t>(numParts, targetNumParts,
1374 numCriteria, partSizes.view(0, numCriteria),
1375 metrics[1].getGlobalSum(), wgts,
1378 metrics[1].setMinImbalance(1.0 + min);
1379 metrics[1].setMaxImbalance(1.0 + max);
1380 metrics[1].setAvgImbalance(1.0 + avg);
1382 if (metrics.size() > 2){
1389 for (
int vdim=0; vdim < numCriteria; vdim++){
1393 if (partSizes[vdim].size() > 0)
1394 psizes = partSizes[vdim].getRawPtr();
1396 computeImbalances<scalar_t, part_t>(numParts, targetNumParts, psizes,
1397 metrics[next].getGlobalSum(), wgts, min, max, avg);
1399 metrics[next].setMinImbalance(1.0 + min);
1400 metrics[next].setMaxImbalance(1.0 + max);
1401 metrics[next].setAvgImbalance(1.0 + avg);
1413 template <
typename scalar_t,
typename part_t>
1415 part_t targetNumParts, part_t numParts, part_t numNonemptyParts,
1418 os <<
"NUMBER OF PARTS IS " << numParts;
1419 if (numNonemptyParts < numParts)
1420 os <<
" (" << numNonemptyParts <<
" of which are non-empty)";
1422 if (targetNumParts != numParts)
1423 os <<
"TARGET NUMBER OF PARTS IS " << targetNumParts << std::endl;
1425 std::string unset(
"unset");
1429 for (
int i=0; i < infoList.size(); i++)
1430 if (infoList[i].
getName() != unset)
1439 template <
typename scalar_t,
typename part_t>
1441 part_t targetNumParts, part_t numParts,
1444 os <<
"NUMBER OF PARTS IS " << numParts;
1446 if (targetNumParts != numParts)
1447 os <<
"TARGET NUMBER OF PARTS IS " << targetNumParts << std::endl;
1449 std::string unset(
"unset");
1453 for (
int i=0; i < infoList.size(); i++)
1454 if (infoList[i].
getName() != unset)
1463 template <
typename scalar_t,
typename part_t>
1465 part_t targetNumParts, part_t numParts, part_t numNonemptyParts,
1468 ArrayView<MetricValues<scalar_t> > infoList(&info, 1);
1469 printMetrics( os, targetNumParts, numParts, numNonemptyParts, infoList);
1475 template <
typename scalar_t,
typename part_t>
1477 part_t targetNumParts, part_t numParts,
1480 ArrayView<GraphMetricValues<scalar_t> > infoList(&info, 1);
1486 template <
typename scalar_t>
1490 size_t dim = weights.size();
1496 scalar_t nweight = 0;
1501 for (
size_t i=0; i <dim; i++)
1502 nweight += weights[i];
1507 for (
size_t i=0; i <dim; i++)
1508 nweight += (weights[i] * weights[i]);
1510 nweight = sqrt(nweight);
1516 nweight = weights[0];
1517 for (
size_t i=1; i <dim; i++)
1518 if (weights[i] > nweight)
1519 nweight = weights[i];
1535 template<
typename lno_t,
typename scalar_t>
1543 Array<scalar_t> vec(dim, 1.0);
1544 for (
size_t i=0; i < dim; i++)
metricOffset
Enumerator for offsets into metric data.
scalar_t getMaxImbalance() const
Get the imbalance of the most imbalanced part. This is what we normally call the imbalance of a parti...
void setName(std::string name)
Set or reset the name.
void computeImbalances(part_t numParts, part_t targetNumParts, const scalar_t *psizes, scalar_t sumVals, const scalar_t *vals, scalar_t &min, scalar_t &max, scalar_t &avg)
Compute the imbalance.
fast typical checks for valid arguments
void normedPartWeights(const RCP< const Environment > &env, part_t numberOfParts, const ArrayView< const part_t > &parts, const ArrayView< StridedData< lno_t, scalar_t > > &vwgts, multiCriteriaNorm mcNorm, scalar_t *weights)
Compute the total weight in each part on this process.
void setAvgImbalance(scalar_t x)
Set the average imbalance of all parts.
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
void setGlobalSum(scalar_t x)
Set the global sum.
scalar_t getGlobalMin() const
Get the global minimum across all parts.
scalar_t getGlobalMax() const
Get the global maximum across all parts.
metricOffset
Enumerator for offsets into metric data.
scalar_t getGlobalAvg() const
Get the average of the sum over all parts.
void setGlobalMax(scalar_t x)
Set the global maximum across parts.
void printMetrics(std::ostream &os, part_t targetNumParts, part_t numParts, part_t numNonemptyParts, const ArrayView< MetricValues< scalar_t > > &infoList)
Print out a header and the values for a list of metrics.
Defines the PartitioningSolution class.
void setGlobalSum(scalar_t x)
Set the global sum.
scalar_t getLocalSum() const
Get the sum on the local process.
sub-steps, each method's entry and exit
scalar_t getAvgImbalance() const
Get the average of the part imbalances.
const std::string & getName() const
Get the name of the item measured.
const std::string & getName() const
Get the name of the item measured.
void setNorm(multiCriteriaNorm normVal)
Set or reset the norm.
GraphMetricValues(std::string mname)
Constructor.
void getStridedStats(const ArrayView< scalar_t > &v, int stride, int offset, scalar_t &min, scalar_t &max, scalar_t &sum)
Find min, max and sum of metric values.
A PartitioningSolution is a solution to a partitioning problem.
void setGlobalAvg(scalar_t x)
Set the global average (sum / numParts).
void localBugAssertion(const char *file, int lineNum, const char *msg, bool ok, AssertionLevel level) const
Test for valid library behavior on local process only.
The StridedData class manages lists of weights or coordinates.
void setMinImbalance(scalar_t x)
Set the imbalance of the least imbalanced part.
scalar_t getGlobalSum() const
Get the global sum for all parts.
The user parameters, debug, timing and memory profiling output objects, and error checking methods...
void objectMetrics(const RCP< const Environment > &env, const RCP< const Comm< int > > &comm, multiCriteriaNorm mcNorm, const RCP< const typename Adapter::base_adapter_t > &ia, const RCP< const PartitioningSolution< Adapter > > &solution, bool useDegreeAsWeight, const RCP< const GraphModel< typename Adapter::base_adapter_t > > &graphModel, typename Adapter::part_t &numParts, typename Adapter::part_t &numNonemptyParts, ArrayRCP< MetricValues< typename Adapter::scalar_t > > &metrics)
Compute imbalance metrics for a distribution.
static void printHeader(std::ostream &os)
Print a standard header.
GraphMetricValues()
Constructor.
void setName(std::string name)
Set or reset the name.
MetricValues(std::string mname)
Constructor.
void printLine(std::ostream &os) const
Print a standard line of data that fits under the header.
GraphModel defines the interface required for graph models.
scalar_t getGlobalMax() const
Get the global maximum across all parts.
void setGlobalMin(scalar_t x)
Set the global minimum across parts.
A class containing the metrics for one measurable item.
multiCriteriaNorm
Enumerator used in code for multicriteria norm choice.
A class containing the metrics for one measurable item.
MetricValues()
Constructor.
void setGlobalMax(scalar_t x)
Set the global maximum across parts.
void globalSumsByPart(const RCP< const Environment > &env, const RCP< const Comm< int > > &comm, const ArrayView< const part_t > &part, int vwgtDim, const ArrayView< StridedData< lno_t, scalar_t > > &vwgts, multiCriteriaNorm mcNorm, part_t &numParts, part_t &numNonemptyParts, ArrayRCP< MetricValues< scalar_t > > &metrics, ArrayRCP< scalar_t > &globalSums)
Given the local partitioning, compute the global sums in each part.
void setLocalSum(scalar_t x)
Set the sum on the local process.
multiCriteriaNorm getNorm()
Get the norm.
Defines the GraphModel interface.
#define Z2_THROW_OUTSIDE_ERROR(env)
Throw an error returned from outside the Zoltan2 library.
void globalWeightedCutsByPart(const RCP< const Environment > &env, const RCP< const Comm< int > > &comm, const RCP< const GraphModel< typename Adapter::base_adapter_t > > &graph, const ArrayView< const typename Adapter::part_t > &part, typename Adapter::part_t &numParts, ArrayRCP< GraphMetricValues< typename Adapter::scalar_t > > &metrics, ArrayRCP< typename Adapter::scalar_t > &globalSums)
Given the local partitioning, compute the global weighted cuts in each part.
static void printHeader(std::ostream &os)
Print a standard header.
scalar_t getMinImbalance() const
Get the imbalance of the least imbalanced part.
scalar_t getGlobalSum() const
Get the global sum for all parts.
void setMaxImbalance(scalar_t x)
Set the imbalance of the worst imbalanced part. This is what we normally call the imbalance of a part...
void printLine(std::ostream &os) const
Print a standard line of data that fits under the header.
This file defines the StridedData class.
scalar_t normedWeight(ArrayView< scalar_t > weights, multiCriteriaNorm norm)
Compute the norm of the vector of weights.