Zoltan2
Zoltan2_EvaluatePartition.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 // @HEADER
46 //
47 // ***********************************************************************
48 //
49 // Zoltan2: A package of combinatorial algorithms for scientific computing
50 // Copyright 2012 Sandia Corporation
51 //
52 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
53 // the U.S. Government retains certain rights in this software.
54 //
55 // Redistribution and use in source and binary forms, with or without
56 // modification, are permitted provided that the following conditions are
57 // met:
58 //
59 // 1. Redistributions of source code must retain the above copyright
60 // notice, this list of conditions and the following disclaimer.
61 //
62 // 2. Redistributions in binary form must reproduce the above copyright
63 // notice, this list of conditions and the following disclaimer in the
64 // documentation and/or other materials provided with the distribution.
65 //
66 // 3. Neither the name of the Corporation nor the names of the
67 // contributors may be used to endorse or promote products derived from
68 // this software without specific prior written permission.
69 //
70 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
71 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
72 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
73 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
74 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
75 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
76 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
77 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
78 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
79 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
80 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
81 //
82 // Questions? Contact Karen Devine (kddevin@sandia.gov)
83 // Erik Boman (egboman@sandia.gov)
84 // Siva Rajamanickam (srajama@sandia.gov)
85 //
86 // ***********************************************************************
87 //
88 // @HEADER
89 
94 #ifndef ZOLTAN2_EVALUATEPARTITION_HPP
95 #define ZOLTAN2_EVALUATEPARTITION_HPP
96 
97 #include <Zoltan2_Metric.hpp>
99 
100 namespace Zoltan2{
101 
109 template <typename Adapter>
111 
112 private:
113 
114  typedef typename Adapter::base_adapter_t base_adapter_t;
115  typedef typename Adapter::lno_t lno_t;
116  typedef typename Adapter::part_t part_t;
117  typedef typename Adapter::scalar_t scalar_t;
118 
120 
121  part_t numGlobalParts_; // desired
122  part_t targetGlobalParts_; // actual
123  part_t numNonEmpty_; // of actual
124 
126  typedef ArrayRCP<RCP<base_metric_type> > base_metric_array_type;
127 
128  base_metric_array_type metricsBase_;
129 
130  void sharedConstructor(const Adapter *ia,
131  ParameterList *p,
132  const RCP<const Comm<int> > &problemComm,
133  const PartitioningSolution<Adapter> *soln,
134  const RCP<const GraphModel
135  <typename Adapter::base_adapter_t> > &graphModel);
136 
137 public:
138 
148  const Adapter *ia,
149  ParameterList *p,
150  const PartitioningSolution<Adapter> *soln,
151  const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel=
152  Teuchos::null):
153  numGlobalParts_(0), targetGlobalParts_(0), numNonEmpty_(0), metricsBase_()
154  {
155  RCP<const Comm<int> > problemComm = DefaultComm<int>::getComm();
156  sharedConstructor(ia, p, problemComm, soln, graphModel);
157  }
158 
159 
170  const Adapter *ia,
171  ParameterList *p,
172  const RCP<const Comm<int> > &problemComm,
173  const PartitioningSolution<Adapter> *soln,
174  const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel=
175  Teuchos::null):
176  numGlobalParts_(0), targetGlobalParts_(0), numNonEmpty_(0), metricsBase_()
177  {
178  sharedConstructor(ia, p, problemComm, soln, graphModel);
179  }
180 
181 #ifdef HAVE_ZOLTAN2_MPI
182 
192  const Adapter *ia,
193  ParameterList *p,
194  MPI_Comm comm,
195  const PartitioningSolution<Adapter> *soln,
196  const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel=
197  Teuchos::null):
198  numGlobalParts_(0), targetGlobalParts_(0), numNonEmpty_(0), metricsBase_()
199  {
200  RCP<Teuchos::OpaqueWrapper<MPI_Comm> > wrapper =
201  Teuchos::opaqueWrapper(comm);
202  RCP<const Comm<int> > problemComm =
203  rcp<const Comm<int> >(new Teuchos::MpiComm<int>(wrapper));
204  sharedConstructor(ia, p, problemComm, soln, graphModel);
205  }
206 #endif
207 
210  base_metric_array_type getAllMetrics() const {
211  return metricsBase_;
212  }
213 
216  std::vector<std::string> getMetricTypes() const {
217  std::set<std::string> temporarySet;
218  std::vector<std::string> returnVector; // Michel did this because he wants
219  // set behavior but preserve the
220  // ordering and a vector is
221  // convenient later since we may
222  // also be looking at the vector
223  // set of all possible names from
224  // BaseClassMetrics
225  // other option is to use only types
226  // that are currently being used
227  for( int n = 0; n < metricsBase_.size(); ++n ) {
228  std::string checkName = metricsBase_[n]->getMetricType();
229  if (temporarySet.find(checkName) == temporarySet.end()) {
230  temporarySet.insert(checkName);
231  returnVector.push_back(checkName);
232  }
233  }
234  return returnVector;
235  }
236 
239  ArrayView<RCP<base_metric_type> > getAllMetricsOfType(std::string metricType) const
240  {
241  // find the beginning and the end of the contiguous block
242  // the list is an ArrayRCP and must preserve any ordering
243  int beginIndex = -1;
244  int sizeOfArrayView = 0;
245  for(typename base_metric_array_type::size_type n = 0;
246  n < metricsBase_.size(); ++n) {
247  if( metricsBase_[n]->getMetricType() == metricType ) {
248  if (beginIndex == -1) {
249  beginIndex = int(n);
250  }
251  ++sizeOfArrayView;
252  }
253  }
254  if (sizeOfArrayView == 0) {
255  return ArrayView<RCP<base_metric_type> >(); // empty array view
256  }
257  return metricsBase_.view(beginIndex, sizeOfArrayView);
258  }
259 
262  scalar_t getObjectCountImbalance() const {
263  ArrayView<RCP<base_metric_type> > metrics =
265  if( metrics.size() <= 0 ) {
266  throw std::logic_error("getObjectCountImbalance() was called "
267  "but no metrics data was generated for " +
268  std::string(IMBALANCE_METRICS_TYPE_NAME) + "." );
269  }
270  return metrics[0]->getMetricValue("maximum imbalance");
271  }
272 
279  scalar_t getNormedImbalance() const{
280  ArrayView<RCP<base_metric_type> > metrics =
282 
283  if( metrics.size() <= 0 ) {
284  throw std::logic_error("getNormedImbalance() was called "
285  "but no metrics data was generated for " +
286  std::string(IMBALANCE_METRICS_TYPE_NAME) + "." );
287  }
288  if( metrics.size() <= 1 ) {
289  throw std::logic_error("getNormedImbalance() was called "
290  "but the normed data does not exist." );
291  }
292  return metrics[1]->getMetricValue("maximum imbalance");
293  }
294 
297  scalar_t getWeightImbalance(int weightIndex) const {
298  // In this case we could have
299  // Option 1
300  // object count
301  // Option 2
302  // object count
303  // weight 0
304  // Option 3
305  // object count
306  // normed imbalance
307  // weight 0
308  // weight 1
309 
310  // if we have multiple weights (meaning array size if 2 or greater) than
311  // the weights begin on index 2
312  // if we have one weight 0 (option 2) then the weights begin on index 1
313  ArrayView<RCP<base_metric_type> > metrics =
315 
316  int weight0IndexStartsAtThisArrayIndex = ( metrics.size() > 2 ) ? 2 : 1;
317  int numberOfWeights = metrics.size() - weight0IndexStartsAtThisArrayIndex;
318  int indexInArray = weight0IndexStartsAtThisArrayIndex + weightIndex;
319  if( metrics.size() <= indexInArray ) {
320  throw std::logic_error("getWeightImbalance was called with weight index "+
321  std::to_string(weightIndex) +
322  " but the maximum weight available for " +
323  std::string(IMBALANCE_METRICS_TYPE_NAME) +
324  " is weight " + std::to_string(numberOfWeights-1) +
325  "." );
326  }
327  return metrics[indexInArray]->getMetricValue("maximum imbalance");
328  }
329 
332  scalar_t getMaxEdgeCut() const{
333  ArrayView<RCP<base_metric_type> > graphMetrics =
335  if( graphMetrics.size() < 1 ) {
336  throw std::logic_error("getMaxEdgeCut() was called "
337  "but no metrics data was generated for " +
338  std::string(GRAPH_METRICS_TYPE_NAME) + "." );
339  }
340  return graphMetrics[0]->getMetricValue("global maximum");
341  }
342 
345  scalar_t getMaxWeightEdgeCut(int weightIndex) const{
346  ArrayView<RCP<base_metric_type> > graphMetrics =
348  int indexInArray = weightIndex + 1; // changed this - it used to start at 0
349  if( graphMetrics.size() <= 1 ) {
350  throw std::logic_error("getMaxWeightEdgeCut was called with "
351  "weight index " + std::to_string(weightIndex) +
352  " but no weights were available for " +
353  std::string(GRAPH_METRICS_TYPE_NAME) + "." );
354  }
355  else if( graphMetrics.size() <= indexInArray ) {
356  // the size() - 2 is because weight 0 starts at array element 1
357  // (so if the array size is 2, the maximum specified weight index is
358  // weight 0 ( 2-2 = 0 )
359  throw std::logic_error("getMaxWeightEdgeCut was called with "
360  "weight index " + std::to_string(weightIndex) +
361  " but the maximum weight available for " +
362  std::string(GRAPH_METRICS_TYPE_NAME) +
363  " is weight " +
364  std::to_string(graphMetrics.size() - 2) + "." );
365  }
366  return graphMetrics[indexInArray]->getMetricValue("global maximum");
367  }
368 
371  scalar_t getTotalEdgeCut() const{
372  ArrayView<RCP<base_metric_type> > graphMetrics =
374  if( graphMetrics.size() < 1 ) {
375  throw std::logic_error("getTotalEdgeCut() was called but no metrics "
376  "data was generated for " +
377  std::string(GRAPH_METRICS_TYPE_NAME) + "." );
378  }
379  return graphMetrics[0]->getMetricValue("global sum");
380  }
381 
384  scalar_t getTotalWeightEdgeCut(int weightIndex) const{
385  ArrayView<RCP<base_metric_type> > graphMetrics =
387  int indexInArray = weightIndex + 1; // changed this; it used to start at 0
388  if( graphMetrics.size() <= 1 ) {
389  // the size() - 2 is because weight 0 starts at array element 1 (so if
390  // the array size is 2, the maximum specified weight index is
391  // weight 0 ( 2-2 = 0 )
392  throw std::logic_error("getTotalWeightEdgeCut was called with "
393  "weight index " + std::to_string(weightIndex) +
394  " but no weights were available for " +
395  std::string(GRAPH_METRICS_TYPE_NAME) + "." );
396  }
397  else if( graphMetrics.size() <= indexInArray ) {
398  throw std::logic_error("getTotalWeightEdgeCut was called with "
399  "weight index " + std::to_string(weightIndex) +
400  " but the maximum weight available for " +
401  std::string(GRAPH_METRICS_TYPE_NAME) +
402  " is weight " +
403  std::to_string(graphMetrics.size() - 2) + "." );
404  }
405  return graphMetrics[indexInArray]->getMetricValue("global sum");
406  }
407 
410  void printMetrics(std::ostream &os, bool bIncludeUnusedTypes = true) const {
411  std::vector<std::string> types =
412  (bIncludeUnusedTypes ?
414  getMetricTypes());
415  for( auto metricType : types ) {
416  printMetrics (os, metricType);
417  }
418  }
419 
422  void printMetrics(std::ostream &os, std::string metricType) const {
423  // Might need changing. The issue here is that they each have unique
424  // parameters to pass.
425  ArrayView<RCP<base_metric_type> > metrics = getAllMetricsOfType(metricType);
426  if (metrics.size() != 0) {
427  // this could be a critical decision - do we want a blank table with
428  // headers when the list is empty - for debugging that is probably better
429  // but it's very messy to have lots of empty tables in the logs
430  if( metricType == GRAPH_METRICS_TYPE_NAME ) {
431  Zoltan2::printGraphMetrics<scalar_t, part_t>(os, targetGlobalParts_,
432  numGlobalParts_, metrics);
433  }
434  else if( metricType == IMBALANCE_METRICS_TYPE_NAME ) {
435  Zoltan2::printImbalanceMetrics<scalar_t, part_t>(os, targetGlobalParts_,
436  numGlobalParts_,
437  numNonEmpty_, metrics);
438  }
439  }
440  }
441 };
442 
443 // sharedConstructor
444 template <typename Adapter>
445 void EvaluatePartition<Adapter>::sharedConstructor(
446  const Adapter *ia,
447  ParameterList *p,
448  const RCP<const Comm<int> > &comm,
449  const PartitioningSolution<Adapter> *soln,
450  const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel
451 )
452 {
453  RCP<const Comm<int> > problemComm;
454  if (comm == Teuchos::null) {
455  problemComm = DefaultComm<int>::getComm();//communicator is Teuchos default
456  } else {
457  problemComm = comm;
458  }
459 
460  RCP<Environment> env;
461 
462  try{
463  env = rcp(new Environment(*p, problemComm));
464  }
466 
467  env->debug(DETAILED_STATUS, std::string("Entering EvaluatePartition"));
468  env->timerStart(MACRO_TIMERS, "Computing metrics");
469 
470  // Parts to which objects are assigned.
471  size_t numLocalObjects = ia->getLocalNumIDs();
472  ArrayRCP<const part_t> parts;
473 
474  if (soln) {
475  // User provided a partitioning solution; use it.
476  parts = arcp(soln->getPartListView(), 0, numLocalObjects, false);
477  env->localInputAssertion(__FILE__, __LINE__, "parts not set",
478  ((numLocalObjects == 0) || soln->getPartListView()), BASIC_ASSERTION);
479  } else {
480  // User did not provide a partitioning solution;
481  // Use input adapter's partition.
482  const part_t *tmp = NULL;
483  ia->getPartsView(tmp);
484  if (tmp != NULL)
485  parts = arcp(tmp, 0, numLocalObjects, false);
486  else {
487  // User has not provided input parts in input adapter
488  part_t *procs = new part_t[numLocalObjects];
489  for (size_t i=0;i<numLocalObjects;i++) procs[i]=problemComm->getRank();
490  parts = arcp(procs, 0, numLocalObjects, true);
491  }
492  }
493  ArrayView<const part_t> partArray = parts(0, numLocalObjects);
494 
495  // When we add parameters for which weights to use, we
496  // should check those here. For now we compute metrics
497  // using all weights.
498 
500  const Teuchos::ParameterEntry *pe = p->getEntryPtr("partitioning_objective");
501  if (pe){
502  std::string strChoice = pe->getValue<std::string>(&strChoice);
503  if (strChoice == std::string("multicriteria_minimize_total_weight"))
504  mcnorm = normMinimizeTotalWeight;
505  else if (strChoice == std::string("multicriteria_minimize_maximum_weight"))
506  mcnorm = normMinimizeMaximumWeight;
507  }
508 
509  const RCP<const base_adapter_t> bia =
510  rcp(dynamic_cast<const base_adapter_t *>(ia), false);
511 
512  try{
513  imbalanceMetrics<Adapter>(env, problemComm, mcnorm, ia, soln, partArray,
514  graphModel,
515  numGlobalParts_, numNonEmpty_, metricsBase_);
516  }
518 
519  if (soln)
520  targetGlobalParts_ = soln->getTargetGlobalNumberOfParts();
521  else
522  targetGlobalParts_ = problemComm->getSize();
523 
524  env->timerStop(MACRO_TIMERS, "Computing metrics");
525 
526  BaseAdapterType inputType = bia->adapterType();
527 
528  if (inputType == GraphAdapterType ||
529  inputType == MatrixAdapterType ||
530  inputType == MeshAdapterType){
531  env->timerStart(MACRO_TIMERS, "Computing graph metrics");
532  // When we add parameters for which weights to use, we
533  // should check those here. For now we compute graph metrics
534  // using all weights.
535 
536  std::bitset<NUM_MODEL_FLAGS> modelFlags;
537 
538  // Create a GraphModel based on input data.
539 
540  RCP<const GraphModel<base_adapter_t> > graph;
541  if (graphModel == Teuchos::null)
542  graph = rcp(new GraphModel<base_adapter_t>(bia, env, problemComm,
543  modelFlags));
544  else
545  graph = graphModel;
546 
547  // compute weighted cuts
548  ArrayRCP<scalar_t> globalSums;
549  try {
550  globalWeightedCutsByPart<Adapter>(env,
551  problemComm, graph, partArray,
552  numGlobalParts_, metricsBase_,
553  globalSums);
554  }
556 
557  env->timerStop(MACRO_TIMERS, "Computing graph metrics");
558  }
559  env->debug(DETAILED_STATUS,
560  std::string("Exiting EvaluatePartition"));
561 }
562 
563 } // namespace Zoltan2
564 
565 #endif
Zoltan2::BaseAdapter< userTypes_t > base_adapter_t
Time an algorithm (or other entity) as a whole.
fast typical checks for valid arguments
void printMetrics(std::ostream &os, bool bIncludeUnusedTypes=true) const
Print all the metrics based on the metric object type.
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
#define IMBALANCE_METRICS_TYPE_NAME
scalar_t getNormedImbalance() const
Return the object normed weight imbalance. Normed imbalance is only valid if there is at least 2 elem...
Defines the PartitioningSolution class.
scalar_t getObjectCountImbalance() const
Return the object count imbalance.
scalar_t getTotalEdgeCut() const
getTotalEdgeCut
sub-steps, each method&#39;s entry and exit
scalar_t getMaxWeightEdgeCut(int weightIndex) const
getMaxWeightEdgeCuts weighted for the specified index
scalar_t getTotalWeightEdgeCut(int weightIndex) const
getTotalWeightEdgeCut weighted for the specified index
A PartitioningSolution is a solution to a partitioning problem.
BaseAdapterType
An enum to identify general types of adapters.
EvaluatePartition(const Adapter *ia, ParameterList *p, const RCP< const Comm< int > > &problemComm, const PartitioningSolution< Adapter > *soln, const RCP< const GraphModel< typename Adapter::base_adapter_t > > &graphModel=Teuchos::null)
Constructor where Teuchos communicator is specified.
The StridedData class manages lists of weights or coordinates.
std::vector< std::string > getMetricTypes() const
GraphModel defines the interface required for graph models.
multiCriteriaNorm
Enumerator used in code for multicriteria norm choice.
EvaluatePartition(const Adapter *ia, ParameterList *p, const PartitioningSolution< Adapter > *soln, const RCP< const GraphModel< typename Adapter::base_adapter_t > > &graphModel=Teuchos::null)
Constructor where communicator is Teuchos default.
ArrayView< RCP< base_metric_type > > getAllMetricsOfType(std::string metricType) const
Return the metric list for types matching the given metric type.
A class that computes and returns quality metrics.
void printMetrics(std::ostream &os, std::string metricType) const
Print all metrics of type metricType based on the metric object type.
scalar_t getWeightImbalance(int weightIndex) const
Return the imbalance for the requested weight.
scalar_t getMaxEdgeCut() const
Return the max cut for the requested weight.
#define GRAPH_METRICS_TYPE_NAME
base_metric_array_type getAllMetrics() const
Return the full metric list which will be mixed types.