Zoltan2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Zoltan2_EvaluateOrdering.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_EVALUATEORDERING_HPP
51 #define ZOLTAN2_EVALUATEORDERING_HPP
52 
55 
56 namespace Zoltan2{
57 
61 template <typename Adapter>
62 class EvaluateOrdering : public EvaluateBaseClass<Adapter> {
63 
64 private:
65  typedef typename Adapter::lno_t lno_t;
66  typedef typename Adapter::gno_t gno_t;
67  typedef typename Adapter::offset_t offset_t;
68  typedef typename Adapter::part_t part_t;
69  typedef typename Adapter::scalar_t scalar_t;
70  typedef typename Adapter::base_adapter_t base_adapter_t;
71 
72  // To do - this is only appropriate for the local ordering
73  // Need to decide how to organize these classes
74  // Do we potentially want local + global in this class
75  // Do we want to eliminate the EvaluateLocalOrdering and EvaluateGlobalOrdering
76  // derived classes? Or perhaps make them completely independent of each other
77  lno_t bandwidth;
78  lno_t envelope;
79  lno_t separatorSize;
80 
81  void sharedConstructor(const Adapter *ia,
82  ParameterList *p,
83  const RCP<const Comm<int> > &problemComm,
84  const LocalOrderingSolution<lno_t> *localSoln,
85  const GlobalOrderingSolution<gno_t> *globalSoln);
86 public:
87 
97  const Adapter *ia,
98  ParameterList *p,
99  const LocalOrderingSolution<lno_t> *localSoln,
100  const GlobalOrderingSolution<gno_t> *globalSoln)
101  {
102  Teuchos::RCP<const Comm<int> > problemComm = Tpetra::getDefaultComm();
103  sharedConstructor(ia, p, problemComm, localSoln, globalSoln);
104  }
105 
116  const Adapter *ia,
117  ParameterList *p,
118  const RCP<const Comm<int> > &problemComm,
119  const LocalOrderingSolution<lno_t> *localSoln,
120  const GlobalOrderingSolution<gno_t> *globalSoln)
121  {
122  sharedConstructor(ia, p, problemComm, localSoln, globalSoln);
123  }
124 
125 #ifdef HAVE_ZOLTAN2_MPI
126 
136  const Adapter *ia,
137  ParameterList *p,
138  MPI_Comm comm,
139  const LocalOrderingSolution<lno_t> *localSoln,
140  const GlobalOrderingSolution<gno_t> *globalSoln)
141  {
142  RCP<Teuchos::OpaqueWrapper<MPI_Comm> > wrapper =
143  Teuchos::opaqueWrapper(comm);
144  RCP<const Comm<int> > problemComm =
145  rcp<const Comm<int> >(new Teuchos::MpiComm<int>(wrapper));
146  sharedConstructor(ia, p, problemComm, localSoln, globalSoln);
147  }
148 #endif
149 
150  lno_t getBandwidth() const { return bandwidth; }
151  lno_t getEnvelope() const { return envelope; }
152  lno_t getSeparatorSize() const { return separatorSize; }
153 
157  virtual void printMetrics(std::ostream &os) const {
158 
159  // To Do - complete this formatting
160  os << "Ordering Metrics" << std::endl;
161  os << std::setw(20) << " " << std::setw(11) << "ordered" << std::endl;
162  os << std::setw(20) << "envelope" << std::setw(11) << std::setprecision(4)
163  << envelope << std::endl;
164  os << std::setw(20) << "bandwidth" << std::setw(11) << std::setprecision(4)
165  << bandwidth << std::endl;
166  }
167 
169  const RCP<const Environment> &env,
170  const RCP<const Comm<int> > &comm,
171  const Adapter *ia,
173  {
174  env->debug(DETAILED_STATUS, "Entering orderingMetrics"); // begin
175 
176  typedef StridedData<lno_t, scalar_t> input_t;
177 
178  // get graph
179  std::bitset<NUM_MODEL_FLAGS> modelFlags;
180  RCP<GraphModel<base_adapter_t> > graph;
181  const RCP<const base_adapter_t> bia =
182  rcp(dynamic_cast<const base_adapter_t *>(ia), false);
183  graph = rcp(new GraphModel<base_adapter_t>(bia,env,comm,modelFlags));
184  ArrayView<const gno_t> Ids;
185  ArrayView<input_t> vwgts;
186  ArrayView<const gno_t> edgeIds;
187  ArrayView<const offset_t> offsets;
188  ArrayView<input_t> wgts;
189  ArrayView<input_t> vtx;
190  graph->getEdgeList(edgeIds, offsets, wgts);
191  lno_t numVertex = graph->getVertexList(Ids, vwgts);
192 
193  lno_t * perm = localSoln->getPermutationView();
194 
195  // print as matrix - this was debugging code which can be deleted later
196  #define MDM
197  #ifdef MDM
198  for( int checkRank = 0; checkRank < comm->getSize(); ++checkRank ) {
199  comm->barrier();
200  if( checkRank == comm->getRank() ) {
201  std::cout << "-----------------------------------------" << std::endl;
202  std::cout << "Inspect rank: " << checkRank << std::endl;
203  std::cout << std::endl;
204  if(numVertex < 30) { // don't spam if it's too many...
205  Array<lno_t> oldMatrix(numVertex*numVertex);
206  Array<lno_t> newMatrix(numVertex*numVertex);
207 
208  // print the solution permutation
209  std::cout << std::endl << "perm: ";
210  for(lno_t n = 0; n < numVertex; ++n) {
211  std::cout << " " << perm[n] << " ";
212  }
213 
214  lno_t * iperm = localSoln->getPermutationView(true);
215  std::cout << std::endl << "iperm: ";
216  for(lno_t n = 0; n < numVertex; ++n) {
217  std::cout << " " << iperm[n] << " ";
218  }
219  std::cout << std::endl;
220  // write 1's to old matrix (original form) and new matrix (using solution)
221  for (lno_t y = 0; y < numVertex; y++) {
222  for (offset_t n = offsets[y]; n < offsets[y+1]; ++n) {
223  lno_t x = static_cast<lno_t>(edgeIds[n]); // to resolve
224  if (x < numVertex && y < numVertex) { // to develop - for MPI this may not be local
225  oldMatrix[x + y*numVertex] = 1;
226  newMatrix[perm[x] + perm[y]*numVertex] = 1;
227  }
228  }
229  }
230 
231  // print oldMatrix
232  std::cout << std::endl << "original graph in matrix form:" << std::endl;
233  for(lno_t y = 0; y < numVertex; ++y) {
234  for(lno_t x = 0; x < numVertex; ++x) {
235  std::cout << " " << oldMatrix[x + y*numVertex];
236  }
237  std::cout << std::endl;
238  }
239 
240  // print newMatrix
241  std::cout << std::endl << "reordered graph in matrix form:" << std::endl;
242  for(lno_t y = 0; y < numVertex; ++y) {
243  for(lno_t x = 0; x < numVertex; ++x) {
244  std::cout << " " << newMatrix[x + y*numVertex];
245  }
246  std::cout << std::endl;
247  }
248  std::cout << std::endl;
249  }
250  }
251 
252  comm->barrier();
253  }
254  #endif // Ends temporary logging which can be deleted later
255 
256  // calculate bandwidth and envelope for unsolved and solved case
257  lno_t bw_right = 0;
258  lno_t bw_left = 0;
259  envelope = 0;
260 
261  for (lno_t j = 0; j < numVertex; j++) {
262  lno_t y = Ids[j];
263  for (offset_t n = offsets[j]; n < offsets[j+1]; ++n) {
264  lno_t x = static_cast<lno_t>(edgeIds[n]); // to resolve
265  if(x < numVertex) {
266  lno_t x2 = perm[x];
267  lno_t y2 = perm[y];
268 
269  // solved bandwidth calculation
270  lno_t delta_right = y2 - x2;
271  if (delta_right > bw_right) {
272  bw_right = delta_right;
273  }
274  lno_t delta_left = y2 - x2;
275  if (delta_left > bw_left) {
276  bw_left = delta_left;
277  }
278 
279  // solved envelope calculation
280  if(delta_right > 0) {
281  envelope += delta_right;
282  }
283  if(delta_left > 0) {
284  envelope += delta_left;
285  }
286  envelope += 1; // need to check this - do we count center?
287  }
288  }
289  }
290 
291  bandwidth = (bw_left + bw_right + 1);
292 
293  // TO DO - No implementation yet for this metric
294  separatorSize = 0;
295 
296  env->debug(DETAILED_STATUS, "Exiting orderingMetrics"); // end
297  }
298 };
299 
300 // sharedConstructor
301 template <typename Adapter>
302 void EvaluateOrdering<Adapter>::sharedConstructor(
303  const Adapter *ia,
304  ParameterList *p,
305  const RCP<const Comm<int> > &comm,
306  const LocalOrderingSolution<lno_t> *localSoln,
307  const GlobalOrderingSolution<gno_t> *globalSoln)
308 {
309  RCP<const Comm<int> > problemComm = (comm == Teuchos::null) ?
310  Tpetra::getDefaultComm() : comm;
311 
312  RCP<Environment> env;
313  try{
314  env = rcp(new Environment(*p, problemComm));
315  }
317 
318  env->debug(DETAILED_STATUS, std::string("Entering EvaluateOrdering"));
319  env->timerStart(MACRO_TIMERS, "Computing ordering metrics");
320 
321  try{
322  // May want to move these into the specific derived classes
323  // But it depends on whether we eventually may have both types and perhaps
324  // want to combine the metrics
325  if(localSoln) {
326  localOrderingMetrics(env, problemComm, ia, localSoln);
327  }
328 
329  if(globalSoln) {
330  throw std::logic_error("EvaluateOrdering not set up for global ordering.");
331  }
332  }
334  env->timerStop(MACRO_TIMERS, "Computing ordering metrics");
335 
336  env->debug(DETAILED_STATUS, std::string("Exiting EvaluateOrdering"));
337 }
338 
339 template <typename Adapter>
340 class EvaluateLocalOrdering : public EvaluateOrdering<Adapter> {
341 private:
342  typedef typename Adapter::lno_t lno_t;
343 
344 public:
353  const Adapter *ia,
354  ParameterList *p,
355  const LocalOrderingSolution<lno_t> *localSoln) :
356  EvaluateOrdering<Adapter>(ia, p, localSoln, nullptr) {}
357 
367  const Adapter *ia,
368  ParameterList *p,
369  const RCP<const Comm<int> > &problemComm,
370  const LocalOrderingSolution<lno_t> *localSoln) :
371  EvaluateOrdering<Adapter>(ia, p, problemComm, localSoln, nullptr) {}
372 
373 #ifdef HAVE_ZOLTAN2_MPI
374 
384  const Adapter *ia,
385  ParameterList *p,
386  MPI_Comm comm,
387  const LocalOrderingSolution<lno_t> *localSoln) :
388  EvaluateOrdering<Adapter>(ia, p, comm, localSoln, nullptr) {}
389 #endif
390 };
391 
392 template <typename Adapter>
393 class EvaluateGlobalOrdering : public EvaluateOrdering<Adapter> {
394 private:
395  typedef typename Adapter::gno_t gno_t;
396 
397 public:
406  const Adapter *ia,
407  ParameterList *p,
408  const GlobalOrderingSolution<gno_t> *globalSoln) :
409  EvaluateOrdering<Adapter>(ia, p, nullptr, globalSoln) {}
410 
420  const Adapter *ia,
421  ParameterList *p,
422  const RCP<const Comm<int> > &problemComm,
423  const GlobalOrderingSolution<gno_t> *globalSoln) :
424  EvaluateOrdering<Adapter>(ia, p, problemComm, nullptr, globalSoln) {}
425 
426 #ifdef HAVE_ZOLTAN2_MPI
427 
437  const Adapter *ia,
438  ParameterList *p,
439  MPI_Comm comm,
440  const GlobalOrderingSolution<gno_t> *globalSoln) :
441  EvaluateOrdering<Adapter>(ia, p, comm, nullptr, globalSoln) {}
442 #endif
443 };
444 
445 } // namespace Zoltan2
446 
447 #endif
Zoltan2::BaseAdapter< userTypes_t > base_adapter_t
Time an algorithm (or other entity) as a whole.
EvaluateGlobalOrdering(const Adapter *ia, ParameterList *p, const RCP< const Comm< int > > &problemComm, const GlobalOrderingSolution< gno_t > *globalSoln)
Constructor where Teuchos communicator is specified.
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
EvaluateOrdering(const Adapter *ia, ParameterList *p, const LocalOrderingSolution< lno_t > *localSoln, const GlobalOrderingSolution< gno_t > *globalSoln)
Constructor where communicator is Teuchos default.
EvaluateOrdering(const Adapter *ia, ParameterList *p, const RCP< const Comm< int > > &problemComm, const LocalOrderingSolution< lno_t > *localSoln, const GlobalOrderingSolution< gno_t > *globalSoln)
Constructor where Teuchos communicator is specified.
Defines the OrderingSolution class.
A class that computes and returns quality metrics. base class for the local and global ordering vers...
sub-steps, each method&#39;s entry and exit
EvaluateLocalOrdering(const Adapter *ia, ParameterList *p, const LocalOrderingSolution< lno_t > *localSoln)
Constructor where communicator is Teuchos default.
SparseMatrixAdapter_t::part_t part_t
virtual void printMetrics(std::ostream &os) const
Print all metrics of type metricType based on the metric object type Note that parent class currently...
lno_t * getPermutationView(bool inverse=false) const
Get pointer to permutation. If inverse = true, return inverse permutation. By default, perm[i] is where new index i can be found in the old ordering. When inverse==true, perm[i] is where old index i can be found in the new ordering.
The StridedData class manages lists of weights or coordinates.
void localOrderingMetrics(const RCP< const Environment > &env, const RCP< const Comm< int > > &comm, const Adapter *ia, const LocalOrderingSolution< typename Adapter::lno_t > *localSoln)
Base class for the EvaluatePartition and EvaluateOrdering classes.
GraphModel defines the interface required for graph models.
EvaluateGlobalOrdering(const Adapter *ia, ParameterList *p, const GlobalOrderingSolution< gno_t > *globalSoln)
Constructor where communicator is Teuchos default.
EvaluateLocalOrdering(const Adapter *ia, ParameterList *p, const RCP< const Comm< int > > &problemComm, const LocalOrderingSolution< lno_t > *localSoln)
Constructor where Teuchos communicator is specified.