Zoltan2
Zoltan2_Metric.hpp
Go to the documentation of this file.
1 // @HEADER
2 //
3 // ***********************************************************************
4 //
5 // Zoltan2: A package of combinatorial algorithms for scientific computing
6 // Copyright 2012 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Karen Devine (kddevin@sandia.gov)
39 // Erik Boman (egboman@sandia.gov)
40 // Siva Rajamanickam (srajama@sandia.gov)
41 //
42 // ***********************************************************************
43 //
44 // @HEADER
45 
51 #ifndef ZOLTAN2_METRIC_HPP
52 #define ZOLTAN2_METRIC_HPP
53 
54 #include <Zoltan2_StridedData.hpp>
56 #include <Tpetra_Map.hpp>
57 #include <Tpetra_CrsMatrix.hpp>
58 #include <Zoltan2_GraphModel.hpp>
59 
60 #include <cmath>
61 #include <iomanip>
62 #include <iostream>
63 #include <vector>
64 
65 namespace Zoltan2{
66 
68 // Classes
70 
74 template <typename scalar_t>
75  class MetricValues{
76 
77 private:
78  void resetValues(){
79  scalar_t *tmp = new scalar_t [evalNumMetrics];
80  memset(tmp, 0, sizeof(scalar_t) * evalNumMetrics);
81  values_ = arcp(tmp, 0, evalNumMetrics, true);
82  }
83  ArrayRCP<scalar_t> values_;
84  std::string metricName_;
85  multiCriteriaNorm mcnorm_; // store "actualNorm + 1"
86 
87 public:
88 
106 };
107 
115 static void printHeader(std::ostream &os);
116 
118 void printLine(std::ostream &os) const;
119 
121 MetricValues(std::string mname) :
122  values_(), metricName_(mname), mcnorm_(multiCriteriaNorm(0)) {
123  resetValues();}
124 
127  values_(), metricName_("unset"), mcnorm_(multiCriteriaNorm(0)) {
128  resetValues();}
129 
131 void setName(std::string name) { metricName_ = name;}
132 
134 void setNorm(multiCriteriaNorm normVal) {
135  mcnorm_ = multiCriteriaNorm(normVal+1);}
136 
138 void setLocalSum(scalar_t x) { values_[evalLocalSum] = x;}
139 
141 void setGlobalSum(scalar_t x) { values_[evalGlobalSum] = x;}
142 
144 void setGlobalMin(scalar_t x) { values_[evalGlobalMin] = x;}
145 
147 void setGlobalMax(scalar_t x) { values_[evalGlobalMax] = x;}
148 
150 void setGlobalAvg(scalar_t x) { values_[evalGlobalAvg] = x;}
151 
153 void setMinImbalance(scalar_t x) { values_[evalMinImbalance] = x;}
154 
158 void setMaxImbalance(scalar_t x) { values_[evalMaxImbalance] = x;}
159 
161 void setAvgImbalance(scalar_t x) { values_[evalAvgImbalance] = x;}
162 
164 const std::string &getName() const { return metricName_; }
165 
168 
170 scalar_t getLocalSum() const { return values_[evalLocalSum];}
171 
173 scalar_t getGlobalSum() const { return values_[evalGlobalSum];}
174 
176 scalar_t getGlobalMin() const { return values_[evalGlobalMin];}
177 
179 scalar_t getGlobalMax() const { return values_[evalGlobalMax];}
180 
182 scalar_t getGlobalAvg() const { return values_[evalGlobalAvg];}
183 
185 scalar_t getMinImbalance() const { return values_[evalMinImbalance];}
186 
190 scalar_t getMaxImbalance() const { return values_[evalMaxImbalance];}
191 
193 scalar_t getAvgImbalance() const { return values_[evalAvgImbalance];}
194 
195 }; // end class
196 
200 template <typename scalar_t>
202 
203 private:
204  void resetValues(){
205  scalar_t *tmp = new scalar_t [evalNumMetrics];
206  memset(tmp, 0, sizeof(scalar_t) * evalNumMetrics);
207  values_ = arcp(tmp, 0, evalNumMetrics, true);
208  }
209  ArrayRCP<scalar_t> values_;
210  std::string metricName_;
211 
212 public:
213 
220 };
221 
224 static void printHeader(std::ostream &os);
225 
227 void printLine(std::ostream &os) const;
228 
230 GraphMetricValues(std::string mname) :
231  values_(), metricName_(mname) {
232  resetValues();}
233 
236  values_(), metricName_("unset") {
237  resetValues();}
238 
240 void setName(std::string name) { metricName_ = name;}
241 
243 void setGlobalSum(scalar_t x) { values_[evalGlobalSum] = x;}
244 
246 void setGlobalMax(scalar_t x) { values_[evalGlobalMax] = x;}
247 
249 const std::string &getName() const { return metricName_; }
250 
252 scalar_t getGlobalSum() const { return values_[evalGlobalSum];}
253 
255 scalar_t getGlobalMax() const { return values_[evalGlobalMax];}
256 
257 }; // end class
258 
259 template <typename scalar_t>
260  void MetricValues<scalar_t>::printLine(std::ostream &os) const
261 {
262  std::string label(metricName_);
263  if (mcnorm_ > 0){
264  multiCriteriaNorm realNorm = multiCriteriaNorm(mcnorm_ - 1);
265  std::ostringstream oss;
266  switch (realNorm){
267  case normMinimizeTotalWeight: // 1-norm = Manhattan norm
268  oss << metricName_ << " (1)";
269  break;
270  case normBalanceTotalMaximum: // 2-norm = sqrt of sum of squares
271  oss << metricName_ << " (2)";
272  break;
273  case normMinimizeMaximumWeight: // inf-norm = maximum norm
274  oss << metricName_ << " (inf)";
275  break;
276  default:
277  oss << metricName_ << " (?)";
278  break;
279  }
280 
281  label = oss.str();
282  }
283 
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) << " ";
289  os << std::setw(6) << std::setprecision(4) << values_[evalMinImbalance];
290  os << std::setw(6) << std::setprecision(4) << values_[evalMaxImbalance];
291  os << std::setw(6) << std::setprecision(4) << values_[evalAvgImbalance];
292  os << std::endl;
293 }
294 
295 template <typename scalar_t>
296  void MetricValues<scalar_t>::printHeader(std::ostream &os)
297 {
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";
302  os << std::endl;
303 
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";
308  os << std::endl;
309 }
310 
311 template <typename scalar_t>
312  void GraphMetricValues<scalar_t>::printLine(std::ostream &os) const
313 {
314  std::string label(metricName_);
315 
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];
319  os << std::endl;
320 }
321 
322 template <typename scalar_t>
324 {
325  os << std::setw(20) << " ";
326  os << std::setw(24) << "----------SUM----------";
327  os << std::endl;
328 
329  os << std::setw(20) << " ";
330  os << std::setw(12) << "sum" << std::setw(12) << "max";
331  os << std::endl;
332 }
333 
335 // Namespace methods to compute metric values
337 
347 template <typename scalar_t>
348  void getStridedStats(const ArrayView<scalar_t> &v, int stride,
349  int offset, scalar_t &min, scalar_t &max, scalar_t &sum)
350 {
351  if (v.size() < 1) return;
352  min = max = sum = v[offset];
353 
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];
357  sum += v[i];
358  }
359 }
360 
369 template <typename scalar_t>
370  void getStridedStats(const ArrayView<scalar_t> &v, int stride,
371  int offset, scalar_t &max, scalar_t &sum)
372 {
373  if (v.size() < 1) return;
374  max = sum = v[offset];
375 
376  for (int i=offset+stride; i < v.size(); i += stride){
377  if (v[i] > max) max = v[i];
378  sum += v[i];
379  }
380 }
381 
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,
405  const ArrayView<StridedData<lno_t, scalar_t> > &vwgts,
406  multiCriteriaNorm mcNorm,
407  scalar_t *weights)
408 {
409  env->localInputAssertion(__FILE__, __LINE__, "parts or weights",
410  numberOfParts > 0 && vwgts.size() > 0, BASIC_ASSERTION);
411 
412  int numObjects = parts.size();
413  int vwgtDim = vwgts.size();
414 
415  memset(weights, 0, sizeof(scalar_t) * numberOfParts);
416 
417  if (numObjects == 0)
418  return;
419 
420  if (vwgtDim == 0) {
421  for (lno_t i=0; i < numObjects; i++){
422  weights[parts[i]]++;
423  }
424  }
425  else if (vwgtDim == 1){
426  for (lno_t i=0; i < numObjects; i++){
427  weights[parts[i]] += vwgts[0][i];
428  }
429  }
430  else{
431  switch (mcNorm){
433  for (int wdim=0; wdim < vwgtDim; wdim++){
434  for (lno_t i=0; i < numObjects; i++){
435  weights[parts[i]] += vwgts[wdim][i];
436  }
437  } // next weight index
438  break;
439 
441  for (lno_t i=0; i < numObjects; i++){
442  scalar_t ssum = 0;
443  for (int wdim=0; wdim < vwgtDim; wdim++)
444  ssum += (vwgts[wdim][i] * vwgts[wdim][i]);
445  weights[parts[i]] += sqrt(ssum);
446  }
447  break;
448 
450  for (lno_t i=0; i < numObjects; i++){
451  scalar_t max = 0;
452  for (int wdim=0; wdim < vwgtDim; wdim++)
453  if (vwgts[wdim][i] > max)
454  max = vwgts[wdim][i];
455  weights[parts[i]] += max;
456  }
457  break;
458 
459  default:
460  env->localBugAssertion(__FILE__, __LINE__, "invalid norm", false,
462  break;
463  }
464  }
465 }
466 
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,
513  int vwgtDim,
514  const ArrayView<StridedData<lno_t, scalar_t> > &vwgts,
515  multiCriteriaNorm mcNorm,
516  part_t &numParts,
517  part_t &numNonemptyParts,
518  ArrayRCP<MetricValues<scalar_t> > &metrics,
519  ArrayRCP<scalar_t> &globalSums)
520 {
521  env->debug(DETAILED_STATUS, "Entering globalSumsByPart");
523  // Initialize return values
524 
525  numParts = numNonemptyParts = 0;
526 
527  int numMetrics = 1; // "object count"
528  if (vwgtDim) numMetrics++; // "normed weight" or "weight 1"
529  if (vwgtDim > 1) numMetrics += vwgtDim; // "weight n"
530 
531  typedef MetricValues<scalar_t> mv_t;
532  mv_t *newMetrics = new mv_t [numMetrics];
533  env->localMemoryAssertion(__FILE__, __LINE__, numMetrics, newMetrics);
534  ArrayRCP<mv_t> metricArray(newMetrics, 0, numMetrics, true);
535 
536  metrics = metricArray;
537 
539  // Figure out the global number of parts in use.
540  // Verify number of vertex weights is the same everywhere.
541 
542  lno_t localNumObj = part.size();
543  part_t localNum[2], globalNum[2];
544  localNum[0] = static_cast<part_t>(vwgtDim);
545  localNum[1] = 0;
546 
547  for (lno_t i=0; i < localNumObj; i++)
548  if (part[i] > localNum[1]) localNum[1] = part[i];
549 
550  try{
551  reduceAll<int, part_t>(*comm, Teuchos::REDUCE_MAX, 2,
552  localNum, globalNum);
553  }
555 
556  env->globalBugAssertion(__FILE__,__LINE__,
557  "inconsistent number of vertex weights",
558  globalNum[0] == localNum[0], DEBUG_MODE_ASSERTION, comm);
559 
560  part_t nparts = globalNum[1] + 1;
561 
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);
566 
568  // Calculate the local totals by part.
569 
570  scalar_t *localBuf = new scalar_t [globalSumSize];
571  env->localMemoryAssertion(__FILE__, __LINE__, globalSumSize, localBuf);
572  memset(localBuf, 0, sizeof(scalar_t) * globalSumSize);
573 
574  scalar_t *obj = localBuf; // # of objects
575 
576  for (lno_t i=0; i < localNumObj; i++)
577  obj[part[i]]++;
578 
579  if (numMetrics > 1){
580 
581  scalar_t *wgt = localBuf + nparts; // single normed weight or weight 1
582  try{
583  normedPartWeights<scalar_t, lno_t, part_t>(env, nparts,
584  part, vwgts, mcNorm, wgt);
585  }
587 
588  // This code assumes the solution has the part ordered the
589  // same way as the user input. (Bug 5891 is resolved.)
590  if (vwgtDim > 1){
591  wgt += nparts; // individual weights
592  for (int vdim = 0; vdim < vwgtDim; vdim++){
593  for (lno_t i=0; i < localNumObj; i++)
594  wgt[part[i]] += vwgts[vdim][i];
595  wgt += nparts;
596  }
597  }
598  }
599 
600  // Metric: local sums on process
601 
602  int next = 0;
603 
604  metrics[next].setName("object count");
605  metrics[next].setLocalSum(localNumObj);
606  next++;
607 
608  if (numMetrics > 1){
609  scalar_t *wgt = localBuf + nparts; // single normed weight or weight 1
610  scalar_t total = 0.0;
611 
612  for (int p=0; p < nparts; p++){
613  total += wgt[p];
614  }
615 
616  if (vwgtDim == 1)
617  metrics[next].setName("weight 1");
618  else
619  metrics[next].setName("normed weight");
620 
621  metrics[next].setLocalSum(total);
622  next++;
623 
624  if (vwgtDim > 1){
625  for (int vdim = 0; vdim < vwgtDim; vdim++){
626  wgt += nparts;
627  total = 0.0;
628  for (int p=0; p < nparts; p++){
629  total += wgt[p];
630  }
631 
632  std::ostringstream oss;
633  oss << "weight " << vdim+1;
634 
635  metrics[next].setName(oss.str());
636  metrics[next].setLocalSum(total);
637  next++;
638  }
639  }
640  }
641 
643  // Obtain global totals by part.
644 
645  try{
646  reduceAll<int, scalar_t>(*comm, Teuchos::REDUCE_SUM, globalSumSize,
647  localBuf, sumBuf);
648  }
650 
651  delete [] localBuf;
652 
654  // Global sum, min, max, and average over all parts
655 
656  obj = sumBuf; // # of objects
657  scalar_t min=0, max=0, sum=0;
658  next = 0;
659 
660  ArrayView<scalar_t> objVec(obj, nparts);
661  getStridedStats<scalar_t>(objVec, 1, 0, min, max, sum);
662 
663  metrics[next].setGlobalMin(min);
664  metrics[next].setGlobalMax(max);
665  metrics[next].setGlobalSum(sum);
666  metrics[next].setGlobalAvg(sum / nparts);
667  next++;
668 
669  if (numMetrics > 1){
670  scalar_t *wgt = sumBuf + nparts; // single normed weight or weight 1
671 
672  ArrayView<scalar_t> normedWVec(wgt, nparts);
673  getStridedStats<scalar_t>(normedWVec, 1, 0, min, max, sum);
674 
675  if (vwgtDim > 1)
676  metrics[next].setNorm(multiCriteriaNorm(mcNorm));
677 
678  metrics[next].setGlobalMin(min);
679  metrics[next].setGlobalMax(max);
680  metrics[next].setGlobalSum(sum);
681  metrics[next].setGlobalAvg(sum / nparts);
682  next++;
683 
684  if (vwgtDim > 1){
685  for (int vdim=0; vdim < vwgtDim; vdim++){
686  wgt += nparts; // individual weights
687  ArrayView<scalar_t> fromVec(wgt, nparts);
688  getStridedStats<scalar_t>(fromVec, 1, 0, min, max, sum);
689 
690  metrics[next].setGlobalMin(min);
691  metrics[next].setGlobalMax(max);
692  metrics[next].setGlobalSum(sum);
693  metrics[next].setGlobalAvg(sum / nparts);
694  next++;
695  }
696  }
697  }
698 
700  // How many parts do we actually have.
701 
702  numParts = nparts;
703  obj = sumBuf; // # of objects
704 
705  /*for (part_t p=nparts-1; p > 0; p--){
706  if (obj[p] > 0) break;
707  numParts--;
708  }*/
709 
710  numNonemptyParts = numParts;
711 
712  for (part_t p=0; p < numParts; p++)
713  if (obj[p] == 0) numNonemptyParts--;
714 
715  env->debug(DETAILED_STATUS, "Exiting globalSumsByPart");
716 }
717 
747 template <typename Adapter>
749  const RCP<const Environment> &env,
750  const RCP<const Comm<int> > &comm,
751  const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graph,
752  const ArrayView<const typename Adapter::part_t> &part,
753  typename Adapter::part_t &numParts,
755  ArrayRCP<typename Adapter::scalar_t> &globalSums)
756 {
757  env->debug(DETAILED_STATUS, "Entering globalWeightedCutsByPart");
759  // Initialize return values
760 
761  numParts = 0;
762 
763  int ewgtDim = graph->getNumWeightsPerEdge();
764 
765  int numMetrics = 1; // "cut count" or "weight 1"
766  if (ewgtDim > 1) numMetrics = ewgtDim; // "weight n"
767 
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;
773  typedef StridedData<lno_t, scalar_t> input_t;
774 
775  typedef GraphMetricValues<scalar_t> mv_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;
779  typedef Tpetra::global_size_t GST;
780  const GST INVALID = Teuchos::OrdinalTraits<GST>::invalid ();
781 
782  using Teuchos::as;
783 
784  mv_t *newMetrics = new mv_t [numMetrics];
785  env->localMemoryAssertion(__FILE__,__LINE__,numMetrics,newMetrics);
786  ArrayRCP<mv_t> metricArray(newMetrics, 0, numMetrics, true);
787 
788  metrics = metricArray;
789 
791  // Figure out the global number of parts in use.
792  // Verify number of vertex weights is the same everywhere.
793 
794  lno_t localNumObj = part.size();
795  part_t localNum[2], globalNum[2];
796  localNum[0] = static_cast<part_t>(ewgtDim);
797  localNum[1] = 0;
798 
799  for (lno_t i=0; i < localNumObj; i++)
800  if (part[i] > localNum[1]) localNum[1] = part[i];
801 
802  try{
803  reduceAll<int, part_t>(*comm, Teuchos::REDUCE_MAX, 2,
804  localNum, globalNum);
805  }
807 
808  env->globalBugAssertion(__FILE__,__LINE__,
809  "inconsistent number of edge weights",
810  globalNum[0] == localNum[0], DEBUG_MODE_ASSERTION, comm);
811 
812  part_t nparts = globalNum[1] + 1;
813 
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);
818 
820  // Calculate the local totals by part.
821 
822  scalar_t *localBuf = new scalar_t [globalSumSize];
823  env->localMemoryAssertion(__FILE__,__LINE__,globalSumSize,localBuf);
824  memset(localBuf, 0, sizeof(scalar_t) * globalSumSize);
825 
826  scalar_t *cut = localBuf; // # of cuts
827 
828  ArrayView<const gno_t> Ids;
829  ArrayView<input_t> vwgts;
830  //size_t nv =
831  graph->getVertexList(Ids, vwgts);
832 
833  ArrayView<const gno_t> edgeIds;
834  ArrayView<const lno_t> offsets;
835  ArrayView<input_t> wgts;
836  //size_t numLocalEdges =
837  graph->getEdgeList(edgeIds, offsets, wgts);
838  // **************************************************************************
839  // *************************** BUILD MAP FOR ADJS ***************************
840  // **************************************************************************
841 
842  RCP<const map_type> vertexMapG;
843 
844  // Build a list of the global vertex ids...
845  gno_t min = std::numeric_limits<gno_t>::max();
846  size_t maxcols = 0;
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;
851  }
852 
853  gno_t gmin;
854  Teuchos::reduceAll<int, gno_t>(*comm,Teuchos::REDUCE_MIN,1,&min,&gmin);
855 
856  //Generate Map for vertex
857  vertexMapG = rcp(new map_type(INVALID, Ids, gmin, comm));
858 
859  // **************************************************************************
860  // ************************** BUILD GRAPH FOR ADJS **************************
861  // **************************************************************************
862 
863  RCP<sparse_matrix_type> adjsMatrix;
864 
865  // Construct Tpetra::CrsGraph objects.
866  adjsMatrix = rcp (new sparse_matrix_type (vertexMapG, 0));
867 
868  Array<part_t> justOneA(maxcols, 1);
869 
870  for (lno_t localElement=0; localElement<localNumObj; ++localElement){
871  // Insert all columns for global row Ids[localElement]
872  size_t ncols = offsets[localElement+1] - offsets[localElement];
873  adjsMatrix->insertGlobalValues(Ids[localElement],
874  edgeIds(offsets[localElement], ncols),
875  justOneA(0, ncols));
876  }
877 
878  //Fill-complete adjs Graph
879  adjsMatrix->fillComplete ();
880 
881  // Compute part
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]);
885  }
886 
887  // Postmultiply adjsMatrix by part
888  adjsMatrix->rightScale(*scaleVec);
889  Array<gno_t> Indices;
890  Array<part_t> Values;
891 
892  if (!ewgtDim) {
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);
899 
900  for (size_t j=0; j < NumEntries; j++)
901  if (part[i] != Values[j])
902  cut[part[i]]++;
903  }
904 
905  // This code assumes the solution has the part ordered the
906  // same way as the user input. (Bug 5891 is resolved.)
907  } else {
908  scalar_t *wgt = localBuf; // weight 1
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);
916 
917  for (size_t j=0; j < NumEntries; j++)
918  if (part[i] != Values[j])
919  wgt[part[i]] += wgts[edim][offsets[i] + j];
920  }
921  wgt += nparts; // individual weights
922  }
923  }
924 
926  // Obtain global totals by part.
927 
928  try{
929  reduceAll<int, scalar_t>(*comm, Teuchos::REDUCE_SUM, globalSumSize,
930  localBuf, sumBuf);
931  }
933 
934  delete [] localBuf;
935 
937  // Global max and sum over all parts
938 
939  cut = sumBuf; // # of cuts
940  scalar_t max=0, sum=0;
941  int next = 0;
942 
943  ArrayView<scalar_t> cutVec(cut, nparts);
944  getStridedStats<scalar_t>(cutVec, 1, 0, max, sum);
945 
946  metrics[next].setName("cut count");
947  metrics[next].setGlobalMax(max);
948  metrics[next].setGlobalSum(sum);
949 
950  if (ewgtDim){
951  scalar_t *wgt = sumBuf; // weight 1
952 
953  for (int edim=0; edim < ewgtDim; edim++){
954  ArrayView<scalar_t> fromVec(wgt, nparts);
955  getStridedStats<scalar_t>(fromVec, 1, 0, max, sum);
956 
957  std::ostringstream oss;
958  oss << "weight " << edim+1;
959 
960  metrics[next].setName(oss.str());
961  metrics[next].setGlobalMax(max);
962  metrics[next].setGlobalSum(sum);
963  next++;
964  wgt += nparts; // individual weights
965  }
966  }
967 
968  numParts = nparts;
969 
970  env->debug(DETAILED_STATUS, "Exiting globalWeightedCutsByPart");
971 }
972 
1003 template <typename scalar_t, typename part_t>
1004  void computeImbalances(part_t numParts, part_t targetNumParts,
1005  const scalar_t *psizes, scalar_t sumVals , const scalar_t *vals,
1006  scalar_t &min, scalar_t &max, scalar_t &avg)
1007 {
1008  min = sumVals;
1009  max = avg = 0;
1010 
1011  if (sumVals <= 0 || targetNumParts < 1 || numParts < 1)
1012  return;
1013 
1014  if (targetNumParts==1 || numParts==1){
1015  min = max = avg = 0; // 0 imbalance
1016  return;
1017  }
1018 
1019  if (!psizes){
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;
1025  avg += tmp;
1026  if (tmp > max) max = tmp;
1027  if (tmp < min) min = tmp;
1028  }
1029  part_t emptyParts = targetNumParts - numParts;
1030  if (emptyParts > 0){
1031  if (max < 1.0)
1032  max = 1.0; // target divided by target
1033  avg += emptyParts;
1034  }
1035  }
1036  else{
1037  for (part_t p=0; p < targetNumParts; p++){
1038  if (psizes[p] > 0){
1039  if (p < numParts){
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;
1044  avg += tmp;
1045  if (tmp > max) max = tmp;
1046  if (tmp < min) min = tmp;
1047  }
1048  else{
1049  if (max < 1.0)
1050  max = 1.0; // target divided by target
1051  avg += 1.0;
1052  }
1053  }
1054  }
1055  }
1056 
1057  avg /= targetNumParts;
1058 }
1059 
1093 template <typename scalar_t, typename part_t>
1094  void computeImbalances(part_t numParts, part_t targetNumParts,
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)
1098 {
1099  min = sumVals;
1100  max = avg = 0;
1101 
1102  if (sumVals <= 0 || targetNumParts < 1 || numParts < 1)
1103  return;
1104 
1105  if (targetNumParts==1 && numParts==1){
1106  min = max = avg = 0; // 0 imbalance
1107  return;
1108  }
1109 
1110  bool allUniformParts = true;
1111  for (int i=0; i < numSizes; i++){
1112  if (psizes[i].size() != 0){
1113  allUniformParts = false;
1114  break;
1115  }
1116  }
1117 
1118  if (allUniformParts){
1119  computeImbalances<scalar_t, part_t>(numParts, targetNumParts, NULL,
1120  sumVals, vals, min, max, avg);
1121  return;
1122  }
1123 
1124  double uniformSize = 1.0 / targetNumParts;
1125  std::vector<double> sizeVec(numSizes);
1126  for (int i=0; i < numSizes; i++){
1127  sizeVec[i] = uniformSize;
1128  }
1129 
1130  for (part_t p=0; p < numParts; p++){
1131 
1132  // If we have objects in parts that should have 0 objects,
1133  // we don't compute an imbalance. It means that other
1134  // parts have too few objects, and the imbalance will be
1135  // picked up there.
1136 
1137  if (p >= targetNumParts)
1138  break;
1139 
1140  // Vector of target amounts: T
1141 
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]);
1148  }
1149  targetNorm = sqrt(targetNorm);
1150 
1151  // If part is supposed to be empty, we don't compute an
1152  // imbalance. Same argument as above.
1153 
1154  if (targetNorm > 0){
1155 
1156  // Vector of actual amounts: A
1157 
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]);
1164  }
1165  actualNorm = sqrt(actualNorm);
1166 
1167  // |A - T| / |T|
1168 
1169  scalar_t imbalance = actualNorm / targetNorm;
1170 
1171  if (imbalance < min)
1172  min = imbalance;
1173  else if (imbalance > max)
1174  max = imbalance;
1175  avg += imbalance;
1176  }
1177  }
1178 
1179  part_t numEmptyParts = 0;
1180 
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;
1186 
1187  if (nonEmptyPart){
1188  // The partitioning has no objects for this part, which
1189  // is supposed to be non-empty.
1190  numEmptyParts++;
1191  }
1192  }
1193 
1194  if (numEmptyParts > 0){
1195  avg += numEmptyParts;
1196  if (max < 1.0)
1197  max = 1.0; // target divided by target
1198  }
1199 
1200  avg /= targetNumParts;
1201 }
1202 
1235 template <typename Adapter>
1237  const RCP<const Environment> &env,
1238  const RCP<const Comm<int> > &comm,
1239  multiCriteriaNorm mcNorm,
1240  const RCP<const typename Adapter::base_adapter_t> &ia,
1241  const RCP<const PartitioningSolution<Adapter> > &solution,
1242  bool useDegreeAsWeight,
1243  const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel,
1244  typename Adapter::part_t &numParts,
1245  typename Adapter::part_t &numNonemptyParts,
1246  ArrayRCP<MetricValues<typename Adapter::scalar_t> > &metrics)
1247 {
1248  env->debug(DETAILED_STATUS, "Entering objectMetrics");
1249 
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;
1255  typedef StridedData<lno_t, scalar_t> sdata_t;
1256 
1257  // Local number of objects.
1258 
1259  size_t numLocalObjects = ia->getLocalNumIDs();
1260 
1261  // Parts to which objects are assigned.
1262 
1263  const part_t *parts;
1264  if (solution != Teuchos::null) {
1265  parts = solution->getPartListView();
1266  env->localInputAssertion(__FILE__, __LINE__, "parts not set",
1267  ((numLocalObjects == 0) || parts), BASIC_ASSERTION);
1268  } else {
1269  part_t *procs = new part_t [numLocalObjects];
1270  for (size_t i = 0; i < numLocalObjects; i++) procs[i] = comm->getRank();
1271  parts = procs;
1272  }
1273  ArrayView<const part_t> partArray(parts, numLocalObjects);
1274 
1275  // Weights, if any, for each object.
1276 
1277  int nWeights = ia->getNumWeightsPerID();
1278  int numCriteria = (nWeights > 0 ? nWeights : 1);
1279  Array<sdata_t> weights(numCriteria);
1280 
1281  if (nWeights == 0){
1282  // One set of uniform weights is implied.
1283  // StridedData default constructor creates length 0 strided array.
1284  weights[0] = sdata_t();
1285  }
1286  else{
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;
1293  graph = rcp(new GraphModel<base_adapter_t>(ia, env, comm, modelFlags));
1294  graph->getVertexList(Ids, vwgts);
1295  } else {
1296  graphModel->getVertexList(Ids, vwgts);
1297  }
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];
1302  }
1303  ArrayRCP<const scalar_t> wgtArray(wgt,0,numLocalObjects,false);
1304  weights[i] = sdata_t(wgtArray, 1);
1305  }
1306  } else {
1307  for (int i=0; i < nWeights; i++){
1308  const scalar_t *wgt;
1309  int stride;
1310  ia->getWeightsView(wgt, stride, i);
1311  ArrayRCP<const scalar_t> wgtArray(wgt,0,stride*numLocalObjects,false);
1312  weights[i] = sdata_t(wgtArray, stride);
1313  }
1314  }
1315  }
1316 
1317  // Relative part sizes, if any, assigned to the parts.
1318 
1319  part_t targetNumParts = solution->getTargetGlobalNumberOfParts();
1320  scalar_t *psizes = NULL;
1321 
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);
1329  }
1330  partSizes[dim] = arcp(psizes, 0, targetNumParts, true);
1331  }
1332  }
1333 
1335  // Get number of parts, and the number that are non-empty.
1336  // Get sums per part of objects, individual weights, and normed weight sums.
1337 
1338  ArrayRCP<scalar_t> globalSums;
1339 
1340  try{
1341  globalSumsByPart<scalar_t, lno_t, part_t>(env, comm,
1342  partArray, nWeights, weights.view(0, numCriteria), mcNorm,
1343  numParts, numNonemptyParts, metrics, globalSums);
1344  }
1346 
1348  // Compute imbalances for the object count.
1349  // (Use first index of part sizes.)
1350 
1351  scalar_t *objCount = globalSums.getRawPtr();
1352  scalar_t min, max, avg;
1353  psizes=NULL;
1354 
1355  if (partSizes[0].size() > 0)
1356  psizes = partSizes[0].getRawPtr();
1357 
1358  computeImbalances<scalar_t, part_t>(numParts, targetNumParts, psizes,
1359  metrics[0].getGlobalSum(), objCount,
1360  min, max, avg);
1361 
1362  metrics[0].setMinImbalance(1.0 + min);
1363  metrics[0].setMaxImbalance(1.0 + max);
1364  metrics[0].setAvgImbalance(1.0 + avg);
1365 
1367  // Compute imbalances for the normed weight sum.
1368 
1369  scalar_t *wgts = globalSums.getRawPtr() + numParts;
1370 
1371  if (metrics.size() > 1){
1372 
1373  computeImbalances<scalar_t, part_t>(numParts, targetNumParts,
1374  numCriteria, partSizes.view(0, numCriteria),
1375  metrics[1].getGlobalSum(), wgts,
1376  min, max, avg);
1377 
1378  metrics[1].setMinImbalance(1.0 + min);
1379  metrics[1].setMaxImbalance(1.0 + max);
1380  metrics[1].setAvgImbalance(1.0 + avg);
1381 
1382  if (metrics.size() > 2){
1383 
1385  // Compute imbalances for each individual weight.
1386 
1387  int next = 2;
1388 
1389  for (int vdim=0; vdim < numCriteria; vdim++){
1390  wgts += numParts;
1391  psizes = NULL;
1392 
1393  if (partSizes[vdim].size() > 0)
1394  psizes = partSizes[vdim].getRawPtr();
1395 
1396  computeImbalances<scalar_t, part_t>(numParts, targetNumParts, psizes,
1397  metrics[next].getGlobalSum(), wgts, min, max, avg);
1398 
1399  metrics[next].setMinImbalance(1.0 + min);
1400  metrics[next].setMaxImbalance(1.0 + max);
1401  metrics[next].setAvgImbalance(1.0 + avg);
1402  next++;
1403  }
1404  }
1405 
1406  }
1407  env->debug(DETAILED_STATUS, "Exiting objectMetrics");
1408 }
1409 
1413 template <typename scalar_t, typename part_t>
1414  void printMetrics( std::ostream &os,
1415  part_t targetNumParts, part_t numParts, part_t numNonemptyParts,
1416  const ArrayView<MetricValues<scalar_t> > &infoList)
1417 {
1418  os << "NUMBER OF PARTS IS " << numParts;
1419  if (numNonemptyParts < numParts)
1420  os << " (" << numNonemptyParts << " of which are non-empty)";
1421  os << std::endl;
1422  if (targetNumParts != numParts)
1423  os << "TARGET NUMBER OF PARTS IS " << targetNumParts << std::endl;
1424 
1425  std::string unset("unset");
1426 
1428 
1429  for (int i=0; i < infoList.size(); i++)
1430  if (infoList[i].getName() != unset)
1431  infoList[i].printLine(os);
1432 
1433  os << std::endl;
1434 }
1435 
1439 template <typename scalar_t, typename part_t>
1440  void printMetrics( std::ostream &os,
1441  part_t targetNumParts, part_t numParts,
1442  const ArrayView<GraphMetricValues<scalar_t> > &infoList)
1443 {
1444  os << "NUMBER OF PARTS IS " << numParts;
1445  os << std::endl;
1446  if (targetNumParts != numParts)
1447  os << "TARGET NUMBER OF PARTS IS " << targetNumParts << std::endl;
1448 
1449  std::string unset("unset");
1450 
1452 
1453  for (int i=0; i < infoList.size(); i++)
1454  if (infoList[i].getName() != unset)
1455  infoList[i].printLine(os);
1456 
1457  os << std::endl;
1458 }
1459 
1463 template <typename scalar_t, typename part_t>
1464  void printMetrics( std::ostream &os,
1465  part_t targetNumParts, part_t numParts, part_t numNonemptyParts,
1466  const MetricValues<scalar_t> &info)
1467 {
1468  ArrayView<MetricValues<scalar_t> > infoList(&info, 1);
1469  printMetrics( os, targetNumParts, numParts, numNonemptyParts, infoList);
1470 }
1471 
1475 template <typename scalar_t, typename part_t>
1476  void printMetrics( std::ostream &os,
1477  part_t targetNumParts, part_t numParts,
1478  const GraphMetricValues<scalar_t> &info)
1479 {
1480  ArrayView<GraphMetricValues<scalar_t> > infoList(&info, 1);
1481  printMetrics( os, targetNumParts, numParts, infoList);
1482 }
1483 
1486 template <typename scalar_t>
1487  scalar_t normedWeight(ArrayView <scalar_t> weights,
1488  multiCriteriaNorm norm)
1489 {
1490  size_t dim = weights.size();
1491  if (dim == 0)
1492  return 0.0;
1493  else if (dim == 1)
1494  return weights[0];
1495 
1496  scalar_t nweight = 0;
1497 
1498  switch (norm){
1501  for (size_t i=0; i <dim; i++)
1502  nweight += weights[i];
1503 
1504  break;
1505 
1507  for (size_t i=0; i <dim; i++)
1508  nweight += (weights[i] * weights[i]);
1509 
1510  nweight = sqrt(nweight);
1511 
1512  break;
1513 
1516  nweight = weights[0];
1517  for (size_t i=1; i <dim; i++)
1518  if (weights[i] > nweight)
1519  nweight = weights[i];
1520 
1521  break;
1522 
1523  default:
1524  Environment env; // default environment
1525  env.localBugAssertion(__FILE__, __LINE__, "invalid norm", false,
1526  BASIC_ASSERTION);
1527  break;
1528  }
1529 
1530  return nweight;
1531 }
1532 
1535 template<typename lno_t, typename scalar_t>
1537  lno_t idx, multiCriteriaNorm norm)
1538 {
1539  size_t dim = weights.size();
1540  if (dim < 1)
1541  return 0;
1542 
1543  Array<scalar_t> vec(dim, 1.0);
1544  for (size_t i=0; i < dim; i++)
1545  if (weights[i].size() > 0)
1546  vec[dim] = weights[i][idx];
1547 
1548  return normedWeight(vec, norm);
1549 }
1550 
1551 } //namespace Zoltan2
1552 #endif
#define INVALID(STR)
metricOffset
Enumerator for offsets into metric data.
checks for logic errors
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.
size_t global_size_t
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.
static ArrayRCP< ArrayRCP< zscalar_t > > weights
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&#39;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.
dictionary vals
Definition: xml2dox.py:186
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.
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.