Zoltan2
Loading...
Searching...
No Matches
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
103
104namespace Zoltan2{
105
113template <typename Adapter>
114class EvaluatePartition : public EvaluateBaseClass<Adapter> {
115
116private:
117
118 typedef typename Adapter::base_adapter_t base_adapter_t;
119 typedef typename Adapter::lno_t lno_t;
120 typedef typename Adapter::part_t part_t;
121 typedef typename Adapter::scalar_t scalar_t;
122
124
125 part_t numGlobalParts_; // desired
126 part_t targetGlobalParts_; // actual
127 part_t numNonEmpty_; // of actual
128
130 typedef ArrayRCP<RCP<base_metric_type> > base_metric_array_type;
131 base_metric_array_type metricsBase_;
132
133protected:
134 void sharedConstructor(const Adapter *ia,
135 ParameterList *p,
136 const RCP<const Comm<int> > &problemComm,
138 const RCP<const GraphModel
139 <typename Adapter::base_adapter_t> > &graphModel);
140
141
142
155 const Adapter *ia,
156 ParameterList *p,
157 const RCP<const Comm<int> > &problemComm,
159 bool force_evaluate,
160 const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel=
161 Teuchos::null):
162 numGlobalParts_(0), targetGlobalParts_(0), numNonEmpty_(0), metricsBase_() {
163 if (force_evaluate){
164 sharedConstructor(ia, p, problemComm, soln, graphModel);
165 }
166
167 }
169 const RCP<const Environment> &_env,
170 const RCP<const Comm<int> > &_problemComm,
171 const RCP<const GraphModel<typename Adapter::base_adapter_t> > &_graph,
172 const ArrayView<const typename Adapter::part_t> &_partArray,
173 typename Adapter::part_t &_numGlobalParts,
174 ArrayRCP<RCP<BaseClassMetrics<typename Adapter::scalar_t> > > &_metricsBase,
175 ArrayRCP<typename Adapter::scalar_t> &_globalSums) {
176 globalWeightedByPart <Adapter>(_env, _problemComm, _graph,
177 _partArray, _numGlobalParts, _metricsBase, _globalSums);
178 }
179public:
181
191 const Adapter *ia,
192 ParameterList *p,
194 const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel=
195 Teuchos::null):
196 numGlobalParts_(0), targetGlobalParts_(0), numNonEmpty_(0), metricsBase_()
197 {
198 Teuchos::RCP<const Comm<int> > problemComm = Tpetra::getDefaultComm();
199 sharedConstructor(ia, p, problemComm, soln, graphModel);
200 }
201
202
213 const Adapter *ia,
214 ParameterList *p,
215 const RCP<const Comm<int> > &problemComm,
217 const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel=
218 Teuchos::null):
219 numGlobalParts_(0), targetGlobalParts_(0), numNonEmpty_(0), metricsBase_()
220 {
221 sharedConstructor(ia, p, problemComm, soln, graphModel);
222 }
223
224#ifdef HAVE_ZOLTAN2_MPI
235 const Adapter *ia,
236 ParameterList *p,
237 MPI_Comm comm,
239 const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel=
240 Teuchos::null):
241 numGlobalParts_(0), targetGlobalParts_(0), numNonEmpty_(0), metricsBase_()
242 {
243 RCP<Teuchos::OpaqueWrapper<MPI_Comm> > wrapper =
244 Teuchos::opaqueWrapper(comm);
245 RCP<const Comm<int> > problemComm =
246 rcp<const Comm<int> >(new Teuchos::MpiComm<int>(wrapper));
247 sharedConstructor(ia, p, problemComm, soln, graphModel);
248 }
249#endif
250
254 // TO DO - note that with current status it probably makes more sense to
255 // break up metricsBase_ into imbalanceMetrics_ and graphMetrics_.
256 // So instead of mixing them just keep two arrays and eliminate this function.
257 // That will clean up several places in this class.
258 ArrayView<RCP<base_metric_type>> getAllMetricsOfType(
259 std::string metricType) const {
260 // find the beginning and the end of the contiguous block
261 // the list is an ArrayRCP and must preserve any ordering
262 int beginIndex = -1;
263 int sizeOfArrayView = 0;
264 for(auto n = 0; n < metricsBase_.size(); ++n) {
265 if( metricsBase_[n]->getMetricType() == metricType ) {
266 if (beginIndex == -1) {
267 beginIndex = int(n);
268 }
269 ++sizeOfArrayView;
270 }
271 }
272 if (sizeOfArrayView == 0) {
273 return ArrayView<RCP<base_metric_type> >(); // empty array view
274 }
275 return metricsBase_.view(beginIndex, sizeOfArrayView);
276 }
277
280 scalar_t getObjectCountImbalance() const {
282 if( metrics.size() <= 0 ) {
283 throw std::logic_error("getObjectCountImbalance() was called "
284 "but no metrics data was generated for " +
285 std::string(IMBALANCE_METRICS_TYPE_NAME) + "." );
286 }
287 return metrics[0]->getMetricValue("maximum imbalance");
288 }
289
296 scalar_t getNormedImbalance() const{
298 if( metrics.size() <= 0 ) {
299 throw std::logic_error("getNormedImbalance() was called "
300 "but no metrics data was generated for " +
301 std::string(IMBALANCE_METRICS_TYPE_NAME) + "." );
302 }
303 if( metrics.size() <= 1 ) {
304 throw std::logic_error("getNormedImbalance() was called "
305 "but the normed data does not exist." );
306 }
307 return metrics[1]->getMetricValue("maximum imbalance");
308 }
309
312 scalar_t getWeightImbalance(int weightIndex) const {
313 // In this case we could have
314 // Option 1
315 // object count
316 // Option 2
317 // object count
318 // weight 0
319 // Option 3
320 // object count
321 // normed imbalance
322 // weight 0
323 // weight 1
324
325 // if we have multiple weights (meaning array size if 2 or greater) than
326 // the weights begin on index 2
327 // if we have one weight 0 (option 2) then the weights begin on index 1
329 int weight0IndexStartsAtThisArrayIndex = ( metrics.size() > 2 ) ? 2 : 1;
330 int numberOfWeights = metrics.size() - weight0IndexStartsAtThisArrayIndex;
331 int indexInArray = weight0IndexStartsAtThisArrayIndex + weightIndex;
332 if( metrics.size() <= indexInArray ) {
333 throw std::logic_error("getWeightImbalance was called with weight index "+
334 std::to_string(weightIndex) +
335 " but the maximum weight available for " +
336 std::string(IMBALANCE_METRICS_TYPE_NAME) +
337 " is weight " + std::to_string(numberOfWeights-1) +
338 "." );
339 }
340 return metrics[indexInArray]->getMetricValue("maximum imbalance");
341 }
342
345 scalar_t getMaxEdgeCut() const{
346 auto graphMetrics = getAllMetricsOfType(GRAPH_METRICS_TYPE_NAME);
347 if( graphMetrics.size() < 1 ) {
348 throw std::logic_error("getMaxEdgeCut() was called "
349 "but no metrics data was generated for " +
350 std::string(GRAPH_METRICS_TYPE_NAME) + "." );
351 }
352 return graphMetrics[0]->getMetricValue("global maximum");
353 }
354
357 scalar_t getMaxWeightEdgeCut(int weightIndex) const{
358 auto graphMetrics = getAllMetricsOfType(GRAPH_METRICS_TYPE_NAME);
359 int indexInArray = weightIndex + 1; // changed this - it used to start at 0
360 if( graphMetrics.size() <= 1 ) {
361 throw std::logic_error("getMaxWeightEdgeCut was called with "
362 "weight index " + std::to_string(weightIndex) +
363 " but no weights were available for " +
364 std::string(GRAPH_METRICS_TYPE_NAME) + "." );
365 }
366 else if( graphMetrics.size() <= indexInArray ) {
367 // the size() - 2 is because weight 0 starts at array element 1
368 // (so if the array size is 2, the maximum specified weight index is
369 // weight 0 ( 2-2 = 0 )
370 throw std::logic_error("getMaxWeightEdgeCut was called with "
371 "weight index " + std::to_string(weightIndex) +
372 " but the maximum weight available for " +
373 std::string(GRAPH_METRICS_TYPE_NAME) +
374 " is weight " +
375 std::to_string(graphMetrics.size() - 2) + "." );
376 }
377 return graphMetrics[indexInArray]->getMetricValue("global maximum");
378 }
379
382 scalar_t getTotalEdgeCut() const{
383 auto graphMetrics = getAllMetricsOfType(GRAPH_METRICS_TYPE_NAME);
384 if( graphMetrics.size() < 1 ) {
385 throw std::logic_error("getTotalEdgeCut() was called but no metrics "
386 "data was generated for " +
387 std::string(GRAPH_METRICS_TYPE_NAME) + "." );
388 }
389 return graphMetrics[0]->getMetricValue("global sum");
390 }
391
394 scalar_t getTotalWeightEdgeCut(int weightIndex) const{
395 auto graphMetrics = getAllMetricsOfType(GRAPH_METRICS_TYPE_NAME);
396 int indexInArray = weightIndex + 1; // changed this; it used to start at 0
397 if( graphMetrics.size() <= 1 ) {
398 // the size() - 2 is because weight 0 starts at array element 1 (so if
399 // the array size is 2, the maximum specified weight index is
400 // weight 0 ( 2-2 = 0 )
401 throw std::logic_error("getTotalWeightEdgeCut was called with "
402 "weight index " + std::to_string(weightIndex) +
403 " but no weights were available for " +
404 std::string(GRAPH_METRICS_TYPE_NAME) + "." );
405 }
406 else if( graphMetrics.size() <= indexInArray ) {
407 throw std::logic_error("getTotalWeightEdgeCut was called with "
408 "weight index " + std::to_string(weightIndex) +
409 " but the maximum weight available for " +
410 std::string(GRAPH_METRICS_TYPE_NAME) +
411 " is weight " +
412 std::to_string(graphMetrics.size() - 2) + "." );
413 }
414 return graphMetrics[indexInArray]->getMetricValue("global sum");
415 }
416
419 scalar_t getTotalMessages() const{
420 auto graphMetrics = getAllMetricsOfType(GRAPH_METRICS_TYPE_NAME);
421 if( graphMetrics.size() < 1 ) {
422 throw std::logic_error("getTotalMessages() was called but no metrics "
423 "data was generated for " +
424 std::string(GRAPH_METRICS_TYPE_NAME) + "." );
425 }
426 // TODO: Would be better to avoid hard coding the array access to [1]
427 return graphMetrics[1]->getMetricValue("global sum");
428 }
429
430
433 scalar_t getMaxMessages() const{
434 auto graphMetrics = getAllMetricsOfType(GRAPH_METRICS_TYPE_NAME);
435 if( graphMetrics.size() < 1 ) {
436 throw std::logic_error("getMaxMessages() was called but no metrics "
437 "data was generated for " +
438 std::string(GRAPH_METRICS_TYPE_NAME) + "." );
439 }
440 // TODO: Would be better to avoid hard coding the array access to [1]
441 return graphMetrics[1]->getMetricValue("global maximum");
442 }
443
446 void printMetrics(std::ostream &os) const {
447 // this could be a critical decision - do we want a blank table with
448 // headers when the list is empty - for debugging that is probably better
449 // but it's very messy to have lots of empty tables in the logs
450 ArrayView<RCP<base_metric_type>> graphMetrics =
452 if (graphMetrics.size() != 0) {
453 Zoltan2::printGraphMetrics<scalar_t, part_t>(os, targetGlobalParts_,
454 numGlobalParts_, graphMetrics);
455 }
456
457 ArrayView<RCP<base_metric_type>> imbalanceMetrics =
459 if (imbalanceMetrics.size() != 0) {
460 Zoltan2::printImbalanceMetrics<scalar_t, part_t>(os, targetGlobalParts_,
461 numGlobalParts_, numNonEmpty_, imbalanceMetrics);
462 }
463 }
464};
465
466// sharedConstructor
467template <typename Adapter>
469 const Adapter *ia,
470 ParameterList *p,
471 const RCP<const Comm<int> > &comm,
473 const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel
474)
475{
476 RCP<const Comm<int> > problemComm;
477 if (comm == Teuchos::null) {
478 problemComm = Tpetra::getDefaultComm();
479 } else {
480 problemComm = comm;
481 }
482
483 RCP<Environment> env;
484
485 try{
486 env = rcp(new Environment(*p, problemComm));
487 }
489
490 env->debug(DETAILED_STATUS, std::string("Entering EvaluatePartition"));
491 env->timerStart(MACRO_TIMERS, "Computing metrics");
492
493 // Parts to which objects are assigned.
494 size_t numLocalObjects = ia->getLocalNumIDs();
495 ArrayRCP<const part_t> parts;
496
497 if (soln) {
498 // User provided a partitioning solution; use it.
499 parts = arcp(soln->getPartListView(), 0, numLocalObjects, false);
500 env->localInputAssertion(__FILE__, __LINE__, "parts not set",
501 ((numLocalObjects == 0) || soln->getPartListView()), BASIC_ASSERTION);
502 } else {
503 // User did not provide a partitioning solution;
504 // Use input adapter's partition.
505 const part_t *tmp = NULL;
506 ia->getPartsView(tmp);
507 if (tmp != NULL)
508 parts = arcp(tmp, 0, numLocalObjects, false);
509 else {
510 // User has not provided input parts in input adapter
511 part_t *procs = new part_t[numLocalObjects];
512 for (size_t i=0;i<numLocalObjects;i++) procs[i]=problemComm->getRank();
513 parts = arcp(procs, 0, numLocalObjects, true);
514 }
515 }
516 ArrayView<const part_t> partArray = parts(0, numLocalObjects);
517
518 // When we add parameters for which weights to use, we
519 // should check those here. For now we compute metrics
520 // using all weights.
521
523 const Teuchos::ParameterEntry *pe = p->getEntryPtr("partitioning_objective");
524 if (pe){
525 std::string strChoice = pe->getValue<std::string>(&strChoice);
526 if (strChoice == std::string("multicriteria_minimize_total_weight"))
528 else if (strChoice == std::string("multicriteria_minimize_maximum_weight"))
530 }
531
532 const RCP<const base_adapter_t> bia =
533 rcp(dynamic_cast<const base_adapter_t *>(ia), false);
534
535 try{
536 imbalanceMetrics<Adapter>(env, problemComm, mcnorm, ia, soln, partArray,
537 graphModel,
538 numGlobalParts_, numNonEmpty_, metricsBase_);
539 }
541
542 if (soln)
543 targetGlobalParts_ = soln->getTargetGlobalNumberOfParts();
544 else
545 targetGlobalParts_ = problemComm->getSize();
546
547 env->timerStop(MACRO_TIMERS, "Computing metrics");
548
549 BaseAdapterType inputType = bia->adapterType();
550
551 if (inputType == GraphAdapterType ||
552 inputType == MatrixAdapterType ||
553 inputType == MeshAdapterType){
554 env->timerStart(MACRO_TIMERS, "Computing graph metrics");
555 // When we add parameters for which weights to use, we
556 // should check those here. For now we compute graph metrics
557 // using all weights.
558
559 std::bitset<NUM_MODEL_FLAGS> modelFlags;
560
561 // Create a GraphModel based on input data.
562
563 RCP<const GraphModel<base_adapter_t> > graph = graphModel;
564 if (graphModel == Teuchos::null) {
565 graph = rcp(new GraphModel<base_adapter_t>(bia, env, problemComm,
566 modelFlags));
567 }
568
569 // compute weighted cuts
570 ArrayRCP<scalar_t> globalSums;
571 try {
572 this->calculate_graph_metrics(env, problemComm, graph, partArray,
573 numGlobalParts_, metricsBase_, globalSums);
574 }
576
577 env->timerStop(MACRO_TIMERS, "Computing graph metrics");
578 }
579
580 env->debug(DETAILED_STATUS, std::string("Exiting EvaluatePartition"));
581}
582
583} // namespace Zoltan2
584
585#endif
Base class for the EvaluatePartition and EvaluateOrdering classes.
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
#define GRAPH_METRICS_TYPE_NAME
#define IMBALANCE_METRICS_TYPE_NAME
Defines the PartitioningSolution class.
The user parameters, debug, timing and memory profiling output objects, and error checking methods.
A class that computes and returns quality metrics.
scalar_t getMaxWeightEdgeCut(int weightIndex) const
getMaxWeightEdgeCuts weighted for the specified index
scalar_t getNormedImbalance() const
Return the object normed weight imbalance. Normed imbalance is only valid if there is at least 2 elem...
void printMetrics(std::ostream &os) const
Print all metrics.
EvaluatePartition(const Adapter *ia, ParameterList *p, const RCP< const Comm< int > > &problemComm, const PartitioningSolution< Adapter > *soln, bool force_evaluate, const RCP< const GraphModel< typename Adapter::base_adapter_t > > &graphModel=Teuchos::null)
Constructor where communicator is Teuchos default, and takes another parameter whether to evaluate me...
scalar_t getTotalMessages() const
getTotalMessages
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.
ArrayView< RCP< base_metric_type > > getAllMetricsOfType(std::string metricType) const
Return the metric list for types matching the given metric type.
scalar_t getObjectCountImbalance() const
Return the object count imbalance.
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.
scalar_t getTotalWeightEdgeCut(int weightIndex) const
getTotalWeightEdgeCut weighted for the specified index
scalar_t getTotalEdgeCut() const
getTotalEdgeCut
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.
virtual void calculate_graph_metrics(const RCP< const Environment > &_env, const RCP< const Comm< int > > &_problemComm, const RCP< const GraphModel< typename Adapter::base_adapter_t > > &_graph, const ArrayView< const typename Adapter::part_t > &_partArray, typename Adapter::part_t &_numGlobalParts, ArrayRCP< RCP< BaseClassMetrics< typename Adapter::scalar_t > > > &_metricsBase, ArrayRCP< typename Adapter::scalar_t > &_globalSums)
void sharedConstructor(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)
scalar_t getMaxMessages() const
getMaxMessages
GraphModel defines the interface required for graph models.
A PartitioningSolution is a solution to a partitioning problem.
const part_t * getPartListView() const
Returns the part list corresponding to the global ID list.
size_t getTargetGlobalNumberOfParts() const
Returns the global number of parts desired in the solution.
The StridedData class manages lists of weights or coordinates.
Created by mbenlioglu on Aug 31, 2020.
@ MACRO_TIMERS
Time an algorithm (or other entity) as a whole.
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.
@ DETAILED_STATUS
sub-steps, each method's entry and exit
@ BASIC_ASSERTION
fast typical checks for valid arguments
BaseAdapterType
An enum to identify general types of adapters.
@ GraphAdapterType
graph data
@ MatrixAdapterType
matrix data
@ MeshAdapterType
mesh data
multiCriteriaNorm
Enumerator used in code for multicriteria norm choice.