Zoltan2
Zoltan2_ImbalanceMetricsUtility.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 
49 #ifndef ZOLTAN2_IMBALANCEMETRICSUTILITY_HPP
50 #define ZOLTAN2_IMBALANCEMETRICSUTILITY_HPP
51 
54 
55 namespace Zoltan2{
56 
100 template <typename scalar_t, typename lno_t, typename part_t>
102  const RCP<const Environment> &env,
103  const RCP<const Comm<int> > &comm,
104  const ArrayView<const part_t> &part,
105  int vwgtDim,
106  const ArrayView<StridedData<lno_t, scalar_t> > &vwgts,
107  multiCriteriaNorm mcNorm,
108  part_t targetNumParts,
109  part_t &numExistingParts,
110  part_t &numNonemptyParts,
111  ArrayRCP<RCP<BaseClassMetrics<scalar_t> > > &metrics,
112  ArrayRCP<scalar_t> &globalSums)
113 {
114  env->debug(DETAILED_STATUS, "Entering globalSumsByPart");
116  // Initialize return values
117 
118  numExistingParts = numNonemptyParts = 0;
119 
120  int numMetrics = 1; // "object count"
121  if (vwgtDim) numMetrics++; // "normed weight" or "weight 0"
122  if (vwgtDim > 1) numMetrics += vwgtDim; // "weight n"
123 
124  // add some more metrics to the array
125  typedef ImbalanceMetrics<scalar_t> mv_t;
126  typedef typename ArrayRCP<RCP<BaseClassMetrics<scalar_t> > >::size_type array_size_type;
127  metrics.resize( metrics.size() + numMetrics );
128  for(array_size_type n = metrics.size()-numMetrics; n < metrics.size(); ++n) {
129  mv_t * newMetric = new mv_t; // allocate the new memory
130 
131  // moved this here because we now allocate the polymorphic classes
132  // we should probably reorganize these functions so all data
133  // setup is done on the derived classes
134  // then as a last step we can insert them into the general array of
135  // MetricBase types
136  if (vwgtDim > 1) {
137  newMetric->setNorm(multiCriteriaNorm(mcNorm));
138  }
139 
140  env->localMemoryAssertion(__FILE__,__LINE__,1,newMetric); // check errors
141  metrics[n] = rcp( newMetric ); // create the new members
142  }
143  array_size_type next = metrics.size() - numMetrics; // MDM - this is most
144  // likely temporary to
145  // preserve the format
146  // here - we are now
147  // filling a larger array
148  // so we may not have
149  // started at 0
150 
152  // Figure out the global number of parts in use.
153  // Verify number of vertex weights is the same everywhere.
154 
155  lno_t localNumObj = part.size();
156  part_t localNum[2], globalNum[2];
157  localNum[0] = static_cast<part_t>(vwgtDim);
158  localNum[1] = 0;
159 
160  for (lno_t i=0; i < localNumObj; i++)
161  if (part[i] > localNum[1]) localNum[1] = part[i];
162 
163  try{
164  reduceAll<int, part_t>(*comm, Teuchos::REDUCE_MAX, 2,
165  localNum, globalNum);
166  }
168 
169  env->globalBugAssertion(__FILE__,__LINE__,
170  "inconsistent number of vertex weights",
171  globalNum[0] == localNum[0], DEBUG_MODE_ASSERTION, comm);
172 
173  part_t maxPartPlusOne = globalNum[1] + 1; // Range of possible part IDs:
174  // [0,maxPartPlusOne)
175 
176  part_t globalSumSize = maxPartPlusOne * numMetrics;
177  scalar_t * sumBuf = new scalar_t [globalSumSize];
178  env->localMemoryAssertion(__FILE__, __LINE__, globalSumSize, sumBuf);
179  globalSums = arcp(sumBuf, 0, globalSumSize);
180 
182  // Calculate the local totals by part.
183 
184  scalar_t *localBuf = new scalar_t [globalSumSize];
185  env->localMemoryAssertion(__FILE__, __LINE__, globalSumSize, localBuf);
186  memset(localBuf, 0, sizeof(scalar_t) * globalSumSize);
187 
188  scalar_t *obj = localBuf; // # of objects
189 
190  for (lno_t i=0; i < localNumObj; i++)
191  obj[part[i]]++;
192 
193  if (numMetrics > 1){
194 
195  scalar_t *wgt = localBuf+maxPartPlusOne; // single normed weight or weight 0
196  try{
197  normedPartWeights<scalar_t, lno_t, part_t>(env, maxPartPlusOne,
198  part, vwgts, mcNorm, wgt);
199  }
201 
202  // This code assumes the solution has the part ordered the
203  // same way as the user input. (Bug 5891 is resolved.)
204  if (vwgtDim > 1){
205  wgt += maxPartPlusOne; // individual weights
206  for (int vdim = 0; vdim < vwgtDim; vdim++){
207  for (lno_t i=0; i < localNumObj; i++)
208  wgt[part[i]] += vwgts[vdim][i];
209  wgt += maxPartPlusOne;
210  }
211  }
212  }
213 
214  // Metric: local sums on process
215  metrics[next]->setName("object count");
216  metrics[next]->setMetricValue("local sum", localNumObj);
217 
218  next++;
219 
220  if (numMetrics > 1){
221  scalar_t *wgt = localBuf+maxPartPlusOne; // single normed weight or weight 0
222  scalar_t total = 0.0;
223 
224  for (int p=0; p < maxPartPlusOne; p++){
225  total += wgt[p];
226  }
227 
228  if (vwgtDim == 1)
229  metrics[next]->setName("weight 0");
230  else
231  metrics[next]->setName("normed weight");
232 
233  metrics[next]->setMetricValue("local sum", total);
234 
235  next++;
236 
237  if (vwgtDim > 1){
238  for (int vdim = 0; vdim < vwgtDim; vdim++){
239  wgt += maxPartPlusOne;
240  total = 0.0;
241  for (int p=0; p < maxPartPlusOne; p++){
242  total += wgt[p];
243  }
244 
245  std::ostringstream oss;
246  oss << "weight " << vdim;
247 
248  metrics[next]->setName(oss.str());
249  metrics[next]->setMetricValue("local sum", total);
250 
251  next++;
252  }
253  }
254  }
255 
257  // Obtain global totals by part.
258 
259  try{
260  reduceAll<int, scalar_t>(*comm, Teuchos::REDUCE_SUM, globalSumSize,
261  localBuf, sumBuf);
262  }
264 
265  delete [] localBuf;
266 
268  // Global sum, min, max, and average over all parts
269 
270  obj = sumBuf; // # of objects
271  scalar_t min=0, max=0, sum=0;
272  next = metrics.size() - numMetrics; // MDM - this is most likely temporary
273  // to preserve the format here - we are
274  // now filling a larger array so we may
275  // not have started at 0
276 
277 
278  ArrayView<scalar_t> objVec(obj, maxPartPlusOne);
279  getStridedStats<scalar_t>(objVec, 1, 0, min, max, sum);
280  if (maxPartPlusOne < targetNumParts)
281  min = scalar_t(0); // Some of the target parts are empty
282 
283  metrics[next]->setMetricValue("global minimum", min);
284  metrics[next]->setMetricValue("global maximum", max);
285  metrics[next]->setMetricValue("global sum", sum);
286  next++;
287 
288  if (numMetrics > 1){
289  scalar_t *wgt = sumBuf + maxPartPlusOne; // single normed weight or weight 0
290 
291  ArrayView<scalar_t> normedWVec(wgt, maxPartPlusOne);
292  getStridedStats<scalar_t>(normedWVec, 1, 0, min, max, sum);
293  if (maxPartPlusOne < targetNumParts)
294  min = scalar_t(0); // Some of the target parts are empty
295 
296  metrics[next]->setMetricValue("global minimum", min);
297  metrics[next]->setMetricValue("global maximum", max);
298  metrics[next]->setMetricValue("global sum", sum);
299  next++;
300 
301  if (vwgtDim > 1){
302  for (int vdim=0; vdim < vwgtDim; vdim++){
303  wgt += maxPartPlusOne; // individual weights
304  ArrayView<scalar_t> fromVec(wgt, maxPartPlusOne);
305  getStridedStats<scalar_t>(fromVec, 1, 0, min, max, sum);
306  if (maxPartPlusOne < targetNumParts)
307  min = scalar_t(0); // Some of the target parts are empty
308 
309  metrics[next]->setMetricValue("global minimum", min);
310  metrics[next]->setMetricValue("global maximum", max);
311  metrics[next]->setMetricValue("global sum", sum);
312  next++;
313  }
314  }
315  }
316 
318  // How many parts do we actually have.
319 
320  numExistingParts = maxPartPlusOne;
321  obj = sumBuf; // # of objects
322 
323  /*for (part_t p=nparts-1; p > 0; p--){
324  if (obj[p] > 0) break;
325  numExistingParts--;
326  }*/
327 
328  numNonemptyParts = numExistingParts;
329 
330  for (part_t p=0; p < numExistingParts; p++)
331  if (obj[p] == 0) numNonemptyParts--;
332 
333  env->debug(DETAILED_STATUS, "Exiting globalSumsByPart");
334 }
335 
366 template <typename Adapter>
368  const RCP<const Environment> &env,
369  const RCP<const Comm<int> > &comm,
370  multiCriteriaNorm mcNorm,
371  const Adapter *ia,
372  const PartitioningSolution<Adapter> *solution,
373  const ArrayView<const typename Adapter::part_t> &partArray,
374  const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel,
375  typename Adapter::part_t &numExistingParts,
376  typename Adapter::part_t &numNonemptyParts,
377  ArrayRCP<RCP<BaseClassMetrics<typename Adapter::scalar_t> > > &metrics)
378 {
379  env->debug(DETAILED_STATUS, "Entering objectMetrics");
380 
381  typedef typename Adapter::scalar_t scalar_t;
382  typedef typename Adapter::gno_t gno_t;
383  typedef typename Adapter::lno_t lno_t;
384  typedef typename Adapter::part_t part_t;
385  typedef typename Adapter::base_adapter_t base_adapter_t;
386  typedef StridedData<lno_t, scalar_t> sdata_t;
387 
388  // Local number of objects.
389 
390  size_t numLocalObjects = ia->getLocalNumIDs();
391 
392  // Weights, if any, for each object.
393 
394  int nWeights = ia->getNumWeightsPerID();
395  int numCriteria = (nWeights > 0 ? nWeights : 1);
396  Array<sdata_t> weights(numCriteria);
397 
398  if (nWeights == 0){
399  // One set of uniform weights is implied.
400  // StridedData default constructor creates length 0 strided array.
401  weights[0] = sdata_t();
402  }
403  else{
404  // whether vertex degree is ever used as vertex weight.
405  enum BaseAdapterType adapterType = ia->adapterType();
406  bool useDegreeAsWeight = false;
407  if (adapterType == GraphAdapterType) {
408  useDegreeAsWeight = reinterpret_cast<const GraphAdapter
409  <typename Adapter::user_t, typename Adapter::userCoord_t> *>(ia)->
410  useDegreeAsWeight(0);
411  } else if (adapterType == MatrixAdapterType) {
412  useDegreeAsWeight = reinterpret_cast<const MatrixAdapter
413  <typename Adapter::user_t, typename Adapter::userCoord_t> *>(ia)->
414  useDegreeAsWeight(0);
415  } else if (adapterType == MeshAdapterType) {
416  useDegreeAsWeight =
417  reinterpret_cast<const MeshAdapter<typename Adapter::user_t> *>(ia)->
418  useDegreeAsWeight(0);
419  }
420  if (useDegreeAsWeight) {
421  ArrayView<const gno_t> Ids;
422  ArrayView<sdata_t> vwgts;
423  if (graphModel == Teuchos::null) {
424  std::bitset<NUM_MODEL_FLAGS> modelFlags;
425  RCP<GraphModel<base_adapter_t> > graph;
426  const RCP<const base_adapter_t> bia =
427  rcp(dynamic_cast<const base_adapter_t *>(ia), false);
428  graph = rcp(new GraphModel<base_adapter_t>(bia,env,comm,modelFlags));
429  graph->getVertexList(Ids, vwgts);
430  } else {
431  graphModel->getVertexList(Ids, vwgts);
432  }
433  scalar_t *wgt = new scalar_t[numLocalObjects];
434  for (int i=0; i < nWeights; i++){
435  for (size_t j=0; j < numLocalObjects; j++) {
436  wgt[j] = vwgts[i][j];
437  }
438  ArrayRCP<const scalar_t> wgtArray(wgt,0,numLocalObjects,false);
439  weights[i] = sdata_t(wgtArray, 1);
440  }
441  } else {
442  for (int i=0; i < nWeights; i++){
443  const scalar_t *wgt;
444  int stride;
445  ia->getWeightsView(wgt, stride, i);
446  ArrayRCP<const scalar_t> wgtArray(wgt,0,stride*numLocalObjects,false);
447  weights[i] = sdata_t(wgtArray, stride);
448  }
449  }
450  }
451 
452  // Relative part sizes, if any, assigned to the parts.
453 
454  part_t targetNumParts = comm->getSize();
455  scalar_t *psizes = NULL;
456 
457  ArrayRCP<ArrayRCP<scalar_t> > partSizes(numCriteria);
458  if (solution) {
459  targetNumParts = solution->getTargetGlobalNumberOfParts();
460  for (int dim=0; dim < numCriteria; dim++){
461  if (solution->criteriaHasUniformPartSizes(dim) != true){
462  psizes = new scalar_t [targetNumParts];
463  env->localMemoryAssertion(__FILE__, __LINE__, targetNumParts, psizes);
464  for (part_t i=0; i < targetNumParts; i++){
465  psizes[i] = solution->getCriteriaPartSize(dim, i);
466  }
467  partSizes[dim] = arcp(psizes, 0, targetNumParts, true);
468  }
469  }
470  }
471 
473  // Get number of parts, and the number that are non-empty.
474  // Get sums per part of objects, individual weights, and normed weight sums.
475 
476  ArrayRCP<scalar_t> globalSums;
477 
478  int initialMetricCount = metrics.size();
479  try{
480  globalSumsByPart<scalar_t, lno_t, part_t>(env, comm,
481  partArray, nWeights, weights.view(0, numCriteria), mcNorm,
482  targetNumParts, numExistingParts, numNonemptyParts, metrics, globalSums);
483  }
485 
486  int addedMetricsCount = metrics.size() - initialMetricCount;
487 
489  // Compute imbalances for the object count.
490  // (Use first index of part sizes.)
491 
492  int index = initialMetricCount;
493 
494  scalar_t *objCount = globalSums.getRawPtr();
495  scalar_t min, max, avg;
496  psizes=NULL;
497 
498  if (partSizes[0].size() > 0)
499  psizes = partSizes[0].getRawPtr();
500 
501  scalar_t gsum = metrics[index]->getMetricValue("global sum");
502  computeImbalances<scalar_t, part_t>(numExistingParts, targetNumParts, psizes,
503  gsum, objCount, min, max, avg);
504 
505  metrics[index]->setMetricValue("global average", gsum / targetNumParts);
506 
507  metrics[index]->setMetricValue("maximum imbalance", 1.0 + max);
508  metrics[index]->setMetricValue("average imbalance", avg);
509 
511  // Compute imbalances for the normed weight sum.
512 
513  scalar_t *wgts = globalSums.getRawPtr() + numExistingParts;
514 
515  if (addedMetricsCount > 1){
516  ++index;
517  gsum = metrics[index]->getMetricValue("global sum");
518  computeImbalances<scalar_t, part_t>(numExistingParts, targetNumParts,
519  numCriteria, partSizes.view(0, numCriteria), gsum, wgts, min, max, avg);
520 
521  metrics[index]->setMetricValue("global average", gsum / targetNumParts);
522 
523  metrics[index]->setMetricValue("maximum imbalance", 1.0 + max);
524  metrics[index]->setMetricValue("average imbalance", avg);
525 
526  if (addedMetricsCount > 2){
527 
529  // Compute imbalances for each individual weight.
530 
531  ++index;
532 
533  for (int vdim=0; vdim < numCriteria; vdim++){
534  wgts += numExistingParts;
535  psizes = NULL;
536 
537  if (partSizes[vdim].size() > 0)
538  psizes = partSizes[vdim].getRawPtr();
539 
540  gsum = metrics[index]->getMetricValue("global sum");
541  computeImbalances<scalar_t, part_t>(numExistingParts, targetNumParts,
542  psizes, gsum, wgts, min, max, avg);
543 
544  metrics[index]->setMetricValue("global average", gsum / targetNumParts);
545 
546  metrics[index]->setMetricValue("maximum imbalance", 1.0 + max);
547  metrics[index]->setMetricValue("average imbalance", avg);
548  index++;
549  }
550  }
551 
552  }
553  env->debug(DETAILED_STATUS, "Exiting objectMetrics");
554 }
555 
558 template <typename scalar_t, typename part_t>
560  std::ostream &os,
561  part_t targetNumParts,
562  part_t numExistingParts,
563  part_t numNonemptyParts)
564 {
565  os << "Imbalance Metrics: (" << numExistingParts << " existing parts)";
566  if (numNonemptyParts < numExistingParts) {
567  os << " (" << numNonemptyParts << " of which are non-empty)";
568  }
569  os << std::endl;
570  if (targetNumParts != numExistingParts) {
571  os << "Target number of parts is " << targetNumParts << std::endl;
572  }
574 }
575 
578 template <typename scalar_t, typename part_t>
580  std::ostream &os,
581  part_t targetNumParts,
582  part_t numExistingParts,
583  part_t numNonemptyParts,
584  const ArrayView<RCP<BaseClassMetrics<scalar_t>>> &infoList)
585 {
586  printImbalanceMetricsHeader<scalar_t, part_t>(os, targetNumParts,
587  numExistingParts,
588  numNonemptyParts);
589  for (int i=0; i < infoList.size(); i++) {
590  if (infoList[i]->getName() != METRICS_UNSET_STRING) {
591  infoList[i]->printLine(os);
592  }
593  }
594  os << std::endl;
595 }
596 
599 template <typename scalar_t, typename part_t>
601  std::ostream &os,
602  part_t targetNumParts,
603  part_t numExistingParts,
604  part_t numNonemptyParts,
605  RCP<BaseClassMetrics<scalar_t>> metricValue)
606 {
607  printImbalanceMetricsHeader<scalar_t, part_t>(os, targetNumParts,
608  numExistingParts,
609  numNonemptyParts);
610  metricValue->printLine(os);
611 }
612 
613 } //namespace Zoltan2
614 
615 
616 #endif
void imbalanceMetrics(const RCP< const Environment > &env, const RCP< const Comm< int > > &comm, multiCriteriaNorm mcNorm, const Adapter *ia, const PartitioningSolution< Adapter > *solution, const ArrayView< const typename Adapter::part_t > &partArray, const RCP< const GraphModel< typename Adapter::base_adapter_t > > &graphModel, typename Adapter::part_t &numExistingParts, typename Adapter::part_t &numNonemptyParts, ArrayRCP< RCP< BaseClassMetrics< typename Adapter::scalar_t > > > &metrics)
Compute imbalance metrics for a distribution.
Zoltan2::BaseAdapter< userTypes_t > base_adapter_t
checks for logic errors
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
MatrixAdapter defines the adapter interface for matrices.
GraphAdapter defines the interface for graph-based user data.
MeshAdapter defines the interface for mesh input.
static ArrayRCP< ArrayRCP< zscalar_t > > weights
sub-steps, each method&#39;s entry and exit
A PartitioningSolution is a solution to a partitioning problem.
static void printHeader(std::ostream &os)
Print a standard header.
scalar_t getCriteriaPartSize(int idx, part_t part) const
Get the size for a given weight index and a given part.
BaseAdapterType
An enum to identify general types of adapters.
The StridedData class manages lists of weights or coordinates.
void setNorm(multiCriteriaNorm normVal)
Set or reset the norm.
void printImbalanceMetrics(std::ostream &os, part_t targetNumParts, part_t numExistingParts, part_t numNonemptyParts, const ArrayView< RCP< BaseClassMetrics< scalar_t >>> &infoList)
Print out list of imbalance metrics.
GraphModel defines the interface required for graph models.
multiCriteriaNorm
Enumerator used in code for multicriteria norm choice.
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 targetNumParts, part_t &numExistingParts, part_t &numNonemptyParts, ArrayRCP< RCP< BaseClassMetrics< scalar_t > > > &metrics, ArrayRCP< scalar_t > &globalSums)
Given the local partitioning, compute the global sums in each part.
bool criteriaHasUniformPartSizes(int idx) const
Determine if balancing criteria has uniform part sizes. (User can specify differing part sizes...
size_t getTargetGlobalNumberOfParts() const
Returns the global number of parts desired in the solution.
void printImbalanceMetricsHeader(std::ostream &os, part_t targetNumParts, part_t numExistingParts, part_t numNonemptyParts)
Print out header info for imbalance metrics.
#define Z2_THROW_OUTSIDE_ERROR(env)
Throw an error returned from outside the Zoltan2 library.
#define METRICS_UNSET_STRING