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 
50 #ifndef ZOLTAN2_EVALUATEPARTITION_HPP
51 #define ZOLTAN2_EVALUATEPARTITION_HPP
52 
53 #include <Zoltan2_Metric.hpp>
55 
56 namespace Zoltan2{
57 
65 template <typename Adapter>
67 
68 private:
69 
70  typedef typename Adapter::lno_t lno_t;
71  typedef typename Adapter::part_t part_t;
72  typedef typename Adapter::scalar_t scalar_t;
74 
75  const RCP<const Environment> env_;
76 
77  part_t numGlobalParts_; // desired
78  part_t targetGlobalParts_; // actual
79  part_t numNonEmpty_; // of actual
80 
81  ArrayRCP<MetricValues<scalar_t> > metrics_;
82  ArrayRCP<const MetricValues<scalar_t> > metricsConst_;
83  ArrayRCP<GraphMetricValues<scalar_t> > graphMetrics_;
84  ArrayRCP<const GraphMetricValues<scalar_t> > graphMetricsConst_;
85 
86 public:
87 
99  EvaluatePartition(const RCP<const Environment> &env,
100  const RCP<const Comm<int> > &problemComm,
101  const RCP<const typename Adapter::base_adapter_t> &ia,
102  const RCP<const PartitioningSolution<Adapter> > &soln,
103  bool useDegreeAsWeight,
104  const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel=
105  Teuchos::null);
106 
117  EvaluatePartition(const RCP<const Environment> &env,
118  const RCP<const Comm<int> > &problemComm,
119  const RCP<const typename Adapter::base_adapter_t> &ia,
120  const RCP<const PartitioningSolution<Adapter> > &soln,
121  const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel=
122  Teuchos::null);
123 
127  ArrayRCP<const MetricValues<scalar_t> > getMetrics() const{
128  //BDD return metricsConst_;
129  if(metricsConst_.is_null()) return metrics_;
130  return metricsConst_;
131  }
132 
136  void getObjectCountImbalance(scalar_t &imbalance) const{
137  imbalance = metrics_[0].getMaxImbalance();
138  }
139 
145  void getNormedImbalance(scalar_t &imbalance) const{
146  if (metrics_.size() > 1)
147  imbalance = metrics_[1].getMaxImbalance();
148  else
149  imbalance = metrics_[0].getMaxImbalance();
150  }
151 
158  void getWeightImbalance(scalar_t &imbalance, int idx=0) const{
159  imbalance = 0;
160  if (metrics_.size() > 2) // idx of multiple weights
161  imbalance = metrics_[idx+2].getMaxImbalance();
162  else if (metrics_.size() == 2) // only one weight
163  imbalance = metrics_[1].getMaxImbalance();
164  else // no weights, return object count imbalance
165  imbalance = metrics_[0].getMaxImbalance();
166  }
167 
170  void printMetrics(std::ostream &os) const {
171  Zoltan2::printMetrics<scalar_t, part_t>(os,
172  targetGlobalParts_, numGlobalParts_, numNonEmpty_,
173  metrics_.view(0, metrics_.size()));
174  }
175 
179  ArrayRCP<const GraphMetricValues<scalar_t> > getGraphMetrics() const{
180  if(graphMetricsConst_.is_null()) return graphMetrics_;
181  return graphMetricsConst_;
182  }
183 
190  void getWeightCut(scalar_t &cut, int idx=0) const{
191  if (graphMetrics_.size() < idx) // idx too high
192  cut = graphMetrics_[graphMetrics_.size()-1].getGlobalMax();
193  else if (idx < 0) // idx too low
194  cut = graphMetrics_[0].getGlobalMax();
195  else // idx weight
196  cut = graphMetrics_[idx].getGlobalMax();
197  }
198 
201  void printGraphMetrics(std::ostream &os) const {
202  Zoltan2::printMetrics<scalar_t, part_t>(os,
203  targetGlobalParts_, numGlobalParts_,
204  graphMetrics_.view(0, graphMetrics_.size()));
205  }
206 };
207 
208 template <typename Adapter>
210  const RCP<const Environment> &env,
211  const RCP<const Comm<int> > &problemComm,
212  const RCP<const typename Adapter::base_adapter_t> &ia,
213  const RCP<const PartitioningSolution<Adapter> > &soln,
214  bool useDegreeAsWeight,
215  const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel):
216  env_(env), numGlobalParts_(0), targetGlobalParts_(0), numNonEmpty_(0),
217  metrics_(), metricsConst_()
218 {
219 
220  env->debug(DETAILED_STATUS, std::string("Entering EvaluatePartition"));
221  env->timerStart(MACRO_TIMERS, "Computing metrics");
222 
223  // When we add parameters for which weights to use, we
224  // should check those here. For now we compute metrics
225  // using all weights.
226 
227  const Teuchos::ParameterList &pl = env->getParameters();
229 
230  const Teuchos::ParameterEntry *pe = pl.getEntryPtr("partitioning_objective");
231 
232  if (pe){
233  std::string strChoice = pe->getValue<std::string>(&strChoice);
234  if (strChoice == std::string("multicriteria_minimize_total_weight"))
235  mcnorm = normMinimizeTotalWeight;
236  else if (strChoice == std::string("multicriteria_minimize_maximum_weight"))
237  mcnorm = normMinimizeMaximumWeight;
238  }
239 
240  try{
241  objectMetrics<Adapter>(env, problemComm, mcnorm, ia, soln,
242  useDegreeAsWeight, graphModel, numGlobalParts_,
243  numNonEmpty_,metrics_);
244  }
246 
247  targetGlobalParts_ = soln->getTargetGlobalNumberOfParts();
248 
249  env->timerStop(MACRO_TIMERS, "Computing metrics");
250  env->debug(DETAILED_STATUS, std::string("Exiting EvaluatePartition"));
251 }
252 
253 template <typename Adapter>
255  const RCP<const Environment> &env,
256  const RCP<const Comm<int> > &problemComm,
257  const RCP<const typename Adapter::base_adapter_t> &ia,
258  const RCP<const PartitioningSolution<Adapter> > &soln,
259  const RCP<const GraphModel<typename Adapter::base_adapter_t> > &graphModel):
260  env_(env), numGlobalParts_(0), targetGlobalParts_(0),
261  graphMetrics_(), graphMetricsConst_()
262 {
263 
264  env->debug(DETAILED_STATUS,
265  std::string("Entering EvaluatePartition"));
266  env->timerStart(MACRO_TIMERS, "Computing graph metrics");
267  // When we add parameters for which weights to use, we
268  // should check those here. For now we compute graph metrics
269  // using all weights.
270 
271  typedef typename Adapter::part_t part_t;
272  typedef typename Adapter::base_adapter_t base_adapter_t;
273 
274  std::bitset<NUM_MODEL_FLAGS> modelFlags;
275 
276  // Create a GraphModel based on input data.
277 
278  RCP<GraphModel<base_adapter_t> > graph;
279  if (graphModel == Teuchos::null)
280  graph = rcp(new GraphModel<base_adapter_t>(ia,env,problemComm,modelFlags));
281 
282  // Local number of objects.
283 
284  size_t numLocalObjects = ia->getLocalNumIDs();
285 
286  // Parts to which objects are assigned.
287 
288  const part_t *parts;
289  if (soln != Teuchos::null) {
290  parts = soln->getPartListView();
291  env->localInputAssertion(__FILE__, __LINE__, "parts not set",
292  ((numLocalObjects == 0) || parts), BASIC_ASSERTION);
293  } else {
294  part_t *procs = new part_t [numLocalObjects];
295  for (size_t i=0; i<numLocalObjects; i++) procs[i] = problemComm->getRank();
296  parts = procs;
297  }
298  ArrayView<const part_t> partArray(parts, numLocalObjects);
299 
300  ArrayRCP<scalar_t> globalSums;
301 
302  if (graphModel == Teuchos::null) {
303  try{
304  globalWeightedCutsByPart<Adapter>(env,
305  problemComm, graph, partArray,
306  numGlobalParts_, graphMetrics_,
307  globalSums);
308  }
310  } else {
311  try{
312  globalWeightedCutsByPart<Adapter>(env,
313  problemComm, graphModel, partArray,
314  numGlobalParts_, graphMetrics_,
315  globalSums);
316  }
318  }
319 
320  targetGlobalParts_ = soln->getTargetGlobalNumberOfParts();
321 
322  env->timerStop(MACRO_TIMERS, "Computing graph metrics");
323  env->debug(DETAILED_STATUS,
324  std::string("Exiting EvaluatePartition"));
325 }
326 
327 } // namespace Zoltan2
328 
329 #endif
Time an algorithm (or other entity) as a whole.
void printMetrics(std::ostream &os) const
Print all the metrics.
fast typical checks for valid arguments
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
Metric class and namespace methods to compute quality metrics.
ArrayRCP< const MetricValues< scalar_t > > getMetrics() const
Return the metric values.
Defines the PartitioningSolution class.
void getObjectCountImbalance(scalar_t &imbalance) const
Return the object count imbalance.
void getNormedImbalance(scalar_t &imbalance) const
Return the object normed weight imbalance.
sub-steps, each method&#39;s entry and exit
EvaluatePartition(const RCP< const Environment > &env, const RCP< const Comm< int > > &problemComm, const RCP< const typename Adapter::base_adapter_t > &ia, const RCP< const PartitioningSolution< Adapter > > &soln, bool useDegreeAsWeight, const RCP< const GraphModel< typename Adapter::base_adapter_t > > &graphModel=Teuchos::null)
Constructor.
A PartitioningSolution is a solution to a partitioning problem.
void printGraphMetrics(std::ostream &os) const
Print all the metrics.
The StridedData class manages lists of weights or coordinates.
ArrayRCP< const GraphMetricValues< scalar_t > > getGraphMetrics() const
Return the graph metric values.
GraphModel defines the interface required for graph models.
multiCriteriaNorm
Enumerator used in code for multicriteria norm choice.
void getWeightImbalance(scalar_t &imbalance, int idx=0) const
Return the imbalance for the requested weight.
void getWeightCut(scalar_t &cut, int idx=0) const
Return the max cut for the requested weight.
A class that computes and returns quality metrics.