Zoltan2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
APFMeshAdapterTest.cpp
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 
53 /**************************************************************/
54 /* Includes */
55 /**************************************************************/
56 
58 #include <Zoltan2_Environment.hpp>
62 
63 // Teuchos includes
64 #include "Teuchos_RCP.hpp"
65 #include "Teuchos_XMLParameterListHelpers.hpp"
66 #include "Teuchos_Hashtable.hpp"
67 
68 // SCOREC includes
69 #ifdef HAVE_ZOLTAN2_PARMA
70 #include <parma.h>
71 #include <apf.h>
72 #include <apfMesh.h>
73 #include <apfMDS.h>
74 #include <apfMesh2.h>
75 #include <PCU.h>
76 #include <gmi_mesh.h>
77 #endif
78 
79 #include <set>
80 
81 using Teuchos::ParameterList;
82 using Teuchos::RCP;
83 using Teuchos::ArrayView;
84 
85 // Computes and prints ghost metrics (takes in a Hyper graph model)
86 template <typename Adapter>
88  typedef typename Adapter::gno_t gno_t;
89  typedef typename Adapter::lno_t lno_t;
90  typedef typename Adapter::scalar_t scalar_t;
91  typedef typename Adapter::offset_t offset_t;
93 
94  ArrayView<const gno_t> Ids;
95  ArrayView<input_t> wgts;
96  mdl.getEdgeList(Ids,wgts);
97  ArrayView<bool> isOwner;
98  mdl.getOwnedList(isOwner);
99  size_t numOwned = mdl.getLocalNumOwnedVertices();
100  ArrayView<const gno_t> pins;
101  ArrayView<const offset_t> offsets;
102  mdl.getPinList(pins,offsets,wgts);
103 
104  std::set<gno_t> gids;
105  for (size_t i=0;i<mdl.getLocalNumVertices();i++) {
106  if (isOwner[i])
107  gids.insert(Ids[i]);
108  }
109  std::set<gno_t> ghosts;
110  gno_t num_ghosts=0;
111  for (size_t i=0;i<mdl.getLocalNumPins();i++) {
112  gno_t pin = pins[i];
113  if (gids.find(pin)==gids.end()) {
114  num_ghosts++;
115  if (ghosts.find(pin)==ghosts.end())
116  ghosts.insert(pin);
117  }
118  }
119  std::cout<< "[METRIC] " << PCU_Comm_Self() << " Total number of ghosts in the hypergraph: " << num_ghosts << "\n"
120  << "[METRIC] " << PCU_Comm_Self() << " Number of unique ghosts: " << ghosts.size() << "\n";
121  gno_t unique_ghosts =ghosts.size();
122  gno_t owned_and_ghosts =unique_ghosts+numOwned;
123  gno_t max_o_and_g,min_o_and_g;
124  gno_t max_ghosts,max_u_ghosts;
125  gno_t min_ghosts,min_u_ghosts;
126  max_ghosts = min_ghosts = num_ghosts;
127  max_u_ghosts = min_u_ghosts = unique_ghosts;
128  max_o_and_g = min_o_and_g = owned_and_ghosts;
129  double avg_ghosts,avg_u_ghosts,avg_o_and_g;
130  PCU_Add_Ints(&num_ghosts,1);
131  PCU_Add_Ints(&unique_ghosts,1);
132  PCU_Add_Ints(&owned_and_ghosts,1);
133  PCU_Max_Ints(&max_ghosts,1);
134  PCU_Max_Ints(&max_u_ghosts,1);
135  PCU_Max_Ints(&max_o_and_g,1);
136  PCU_Min_Ints(&min_ghosts,1);
137  PCU_Min_Ints(&min_u_ghosts,1);
138  PCU_Min_Ints(&min_o_and_g,1);
139  avg_ghosts = num_ghosts*1.0/PCU_Comm_Peers();
140  avg_u_ghosts = unique_ghosts*1.0/PCU_Comm_Peers();
141  avg_o_and_g = owned_and_ghosts*1.0/PCU_Comm_Peers();
142  if (!PCU_Comm_Self())
143  std::cout<< "[METRIC] Global ghosts in the hypergraph (tot max min avg imb): "
144  << num_ghosts<<" "<<max_ghosts<<" "<<min_ghosts<<" "<<avg_ghosts<<" "
145  <<max_ghosts/avg_ghosts << "\n"
146  << "[METRIC] Global unique ghosts (tot max min avg imb): "
147  << unique_ghosts<<" "<<max_u_ghosts<<" "<<min_u_ghosts<<" "<<avg_u_ghosts<<" "
148  <<max_u_ghosts/avg_u_ghosts << "\n"
149  << "[METRIC] Global owned and ghosts (tot max min avg imb): "
150  << owned_and_ghosts<<" "<<max_o_and_g<<" "<<min_o_and_g<<" "<<avg_o_and_g<<" "
151  <<max_o_and_g/avg_o_and_g << "\n";
152 
153 }
154 
155 /*****************************************************************************/
156 /******************************** MAIN ***************************************/
157 /*****************************************************************************/
158 
159 int main(int narg, char *arg[]) {
160 
161  Tpetra::ScopeGuard tscope(&narg, &arg);
162  Teuchos::RCP<const Teuchos::Comm<int> > CommT = Tpetra::getDefaultComm();
163 
164  int me = CommT->getRank();
165  //int numProcs = CommT->getSize();
166 
167  if (me == 0){
168  std::cout
169  << "====================================================================\n"
170  << "| |\n"
171  << "| Example: Partition APF Mesh |\n"
172  << "| |\n"
173  << "| Questions? Contact Karen Devine (kddevin@sandia.gov), |\n"
174  << "| Erik Boman (egboman@sandia.gov), |\n"
175  << "| Siva Rajamanickam (srajama@sandia.gov). |\n"
176  << "| |\n"
177  << "| Zoltan2's website: http://trilinos.sandia.gov/packages/zoltan2 |\n"
178  << "| Trilinos website: http://trilinos.sandia.gov |\n"
179  << "| |\n"
180  << "====================================================================\n";
181  }
182 
183 
184 #ifdef HAVE_MPI
185  if (me == 0) {
186  std::cout << "PARALLEL executable \n";
187  }
188 #else
189  if (me == 0) {
190  std::cout << "SERIAL executable \n";
191  }
192 #endif
193 
194  /***************************************************************************/
195  /******************************* GET INPUTS ********************************/
196  /***************************************************************************/
197 
198  // default values for command-line arguments
199  std::string meshFileName("4/");
200  std::string modelFileName("torus.dmg");
201  std::string action("parma");
202  std::string parma_method("VtxElm");
203  std::string output_loc("");
204  int nParts = CommT->getSize();
205  double imbalance = 1.1;
206  int layers=2;
207  int ghost_metric=false;
208  // Read run-time options.
209  Teuchos::CommandLineProcessor cmdp (false, false);
210  cmdp.setOption("meshfile", &meshFileName,
211  "Mesh file with APF specifications (.smb file(s))");
212  cmdp.setOption("modelfile", &modelFileName,
213  "Model file with APF specifications (.dmg file)");
214  cmdp.setOption("action", &action,
215  "Method to use: mj, scotch, zoltan_rcb, parma or color");
216  cmdp.setOption("parma_method", &parma_method,
217  "Method to use: Vertex, Element, VtxElm, VtxEdgeElm, Ghost, Shape, or Centroid ");
218  cmdp.setOption("nparts", &nParts,
219  "Number of parts to create");
220  cmdp.setOption("imbalance", &imbalance,
221  "Target imbalance for the partitioning method");
222  cmdp.setOption("output", &output_loc,
223  "Location of new partitioned apf mesh. Ex: 4/torus.smb");
224  cmdp.setOption("layers", &layers,
225  "Number of layers for ghosting");
226  cmdp.setOption("ghost_metric", &ghost_metric,
227  "0 does not compute ghost metric otherwise compute both before and after");
228  cmdp.parse(narg, arg);
229 
230  /***************************************************************************/
231  /********************** GET CELL TOPOLOGY **********************************/
232  /***************************************************************************/
233 
234  // Get dimensions
235  //int dim = 3;
236 
237  /***************************************************************************/
238  /***************************** GENERATE MESH *******************************/
239  /***************************************************************************/
240 
241 #ifdef HAVE_ZOLTAN2_PARMA
242 
243  if (me == 0) std::cout << "Generating mesh ... \n\n";
244 
245  //Setup for SCOREC
246  PCU_Comm_Init();
247 
248  // Generate mesh with MDS
249  gmi_register_mesh();
250  apf::Mesh2* m = apf::loadMdsMesh(modelFileName.c_str(),meshFileName.c_str());
251  apf::verify(m);
252 
253  //Data for APF MeshAdapter
254  std::string primary="region";
255  std::string adjacency="face";
256  if (m->getDimension()==2) {
257  primary="face";
258  adjacency="edge";
259  }
260  bool needSecondAdj=false;
261 
262  // Set parameters for partitioning
263  if (me == 0) std::cout << "Creating parameter list ... \n\n";
264 
265  Teuchos::ParameterList params("test params");
266  params.set("timer_output_stream" , "std::cout");
267 
268  bool do_partitioning = false;
269  if (action == "mj") {
270  do_partitioning = true;
271  params.set("debug_level", "basic_status");
272  params.set("imbalance_tolerance", imbalance);
273  params.set("num_global_parts", nParts);
274  params.set("algorithm", "multijagged");
275  params.set("rectilinear", "yes");
276  }
277  else if (action == "scotch") {
278  do_partitioning = true;
279  params.set("debug_level", "no_status");
280  params.set("imbalance_tolerance", imbalance);
281  params.set("num_global_parts", nParts);
282  params.set("partitioning_approach", "partition");
283  params.set("objects_to_partition","mesh_elements");
284  params.set("algorithm", "scotch");
285  needSecondAdj=true;
286  }
287  else if (action == "zoltan_rcb") {
288  do_partitioning = true;
289  params.set("debug_level", "verbose_detailed_status");
290  params.set("imbalance_tolerance", imbalance);
291  params.set("num_global_parts", nParts);
292  params.set("partitioning_approach", "partition");
293  params.set("algorithm", "zoltan");
294  }
295  else if (action == "parma") {
296  do_partitioning = true;
297  params.set("debug_level", "no_status");
298  params.set("imbalance_tolerance", imbalance);
299  params.set("algorithm", "parma");
300  Teuchos::ParameterList &pparams = params.sublist("parma_parameters",false);
301  pparams.set("parma_method",parma_method);
302  pparams.set("step_size",1.1);
303  if (parma_method=="Ghost") {
304  pparams.set("ghost_layers",layers);
305  pparams.set("ghost_bridge",m->getDimension()-1);
306  }
307  adjacency="vertex";
308  }
309  else if (action=="zoltan_hg") {
310  do_partitioning = true;
311  params.set("debug_level", "no_status");
312  params.set("imbalance_tolerance", imbalance);
313  params.set("algorithm", "zoltan");
314  params.set("num_global_parts", nParts);
315  Teuchos::ParameterList &zparams = params.sublist("zoltan_parameters",false);
316  zparams.set("LB_METHOD","HYPERGRAPH");
317  zparams.set("LB_APPROACH","PARTITION");
318  //params.set("compute_metrics","yes");
319  adjacency="vertex";
320  }
321  else if (action=="hg_ghost") {
322  do_partitioning = true;
323  params.set("debug_level", "no_status");
324  params.set("imbalance_tolerance", imbalance);
325  params.set("algorithm", "zoltan");
326  params.set("num_global_parts", nParts);
327  params.set("hypergraph_model_type","ghosting");
328  params.set("ghost_layers",layers);
329  Teuchos::ParameterList &zparams = params.sublist("zoltan_parameters",false);
330  zparams.set("LB_METHOD","HYPERGRAPH");
331  zparams.set("LB_APPROACH","PARTITION");
332  zparams.set("PHG_EDGE_SIZE_THRESHOLD", "1.0");
333  primary="vertex";
334  adjacency="edge";
335  needSecondAdj=true;
336  }
337  else if (action == "color") {
338  params.set("debug_level", "verbose_detailed_status");
339  params.set("debug_output_file", "kdd");
340  params.set("debug_procs", "all");
341  }
342  Parma_PrintPtnStats(m,"before");
343 
344  // Creating mesh adapter
345  if (me == 0) std::cout << "Creating mesh adapter ... \n\n";
346  typedef Zoltan2::APFMeshAdapter<apf::Mesh2*> inputAdapter_t;
348  typedef Zoltan2::MeshAdapter<apf::Mesh2*> baseMeshAdapter_t;
349 
350  double time_1=PCU_Time();
351  inputAdapter_t *ia =
352  new inputAdapter_t(*CommT, m,primary,adjacency,needSecondAdj);
353  double time_2=PCU_Time();
354 
355 
356  inputAdapter_t::scalar_t* arr =
357  new inputAdapter_t::scalar_t[ia->getLocalNumOf(ia->getPrimaryEntityType())];
358  for (size_t i=0;i<ia->getLocalNumOf(ia->getPrimaryEntityType());i++) {
359  arr[i]=PCU_Comm_Self()+1;
360  }
361 
362  const inputAdapter_t::scalar_t* weights=arr;
363  ia->setWeights(ia->getPrimaryEntityType(),weights,1);
364 
365 
366  if (ghost_metric) {
367  const baseMeshAdapter_t *base_ia = dynamic_cast<const baseMeshAdapter_t*>(ia);
368  Zoltan2::modelFlag_t graphFlags_;
369  RCP<Zoltan2::Environment> env;
370  try{
371  env = rcp(new Zoltan2::Environment(params, Tpetra::getDefaultComm()));
372  }
374 
375  RCP<const Zoltan2::Environment> envConst = Teuchos::rcp_const_cast<const Zoltan2::Environment>(env);
376 
377  RCP<const baseMeshAdapter_t> baseInputAdapter_(base_ia,false);
378  Zoltan2::HyperGraphModel<inputAdapter_t> model(baseInputAdapter_,envConst,CommT,
379  graphFlags_,Zoltan2::HYPEREDGE_CENTRIC);
380  PrintGhostMetrics(model);
381  }
382 
383  // create Partitioning problem
384  double time_3 = PCU_Time();
385  if (do_partitioning) {
386  if (me == 0) std::cout << "Creating partitioning problem ... \n\n";
387 
388  Zoltan2::PartitioningProblem<inputAdapter_t> problem(ia, &params, CommT);
389 
390  // call the partitioner
391  if (me == 0) std::cout << "Calling the partitioner ... \n\n";
392 
393  problem.solve();
394 
395 
396 
397  if (me==0) std::cout << "Applying Solution to Mesh\n\n";
398  apf::Mesh2** new_mesh = &m;
399  ia->applyPartitioningSolution(m,new_mesh,problem.getSolution());
400 
401  // create metric object
402  RCP<quality_t> metricObject =
403  rcp(new quality_t(ia, &params, CommT, &problem.getSolution()));
404 
405  if (!me) {
406  metricObject->printMetrics(std::cout);
407  }
408  }
409  else {
410  if (me == 0) std::cout << "Creating coloring problem ... \n\n";
411 
412  Zoltan2::ColoringProblem<inputAdapter_t> problem(ia, &params);
413 
414  // call the partitioner
415  if (me == 0) std::cout << "Calling the coloring algorithm ... \n\n";
416 
417  problem.solve();
418 
419  problem.printTimers();
420 
421 
422  }
423 
424  double time_4=PCU_Time();
425 
426  //Destroy the adapter
427  ia->destroy();
428  delete [] arr;
429  //Parma_PrintPtnStats(m,"after");
430 
431  if (ghost_metric) {
432  inputAdapter_t ia2(*CommT, m,primary,adjacency,true);
433  const baseMeshAdapter_t *base_ia = dynamic_cast<const baseMeshAdapter_t*>(&ia2);
434 
435  Zoltan2::modelFlag_t graphFlags_;
436  RCP<Zoltan2::Environment> env;
437  try{
438  env = rcp(new Zoltan2::Environment(params, Tpetra::getDefaultComm()));
439  }
441  RCP<const Zoltan2::Environment> envConst = Teuchos::rcp_const_cast<const Zoltan2::Environment>(env);
442  RCP<const baseMeshAdapter_t> baseInputAdapter_(base_ia,false);
443  Zoltan2::HyperGraphModel<inputAdapter_t> model(baseInputAdapter_, envConst, CommT,
444  graphFlags_,Zoltan2::HYPEREDGE_CENTRIC);
445 
446  PrintGhostMetrics(model);
447  ia2.destroy();
448  }
449 
450  if (output_loc!="") {
451  m->writeNative(output_loc.c_str());
452  }
453 
454  // delete mesh
455  if (me == 0) std::cout << "Deleting the mesh ... \n\n";
456  time_4-=time_3;
457  time_2-=time_1;
458  PCU_Max_Doubles(&time_2,1);
459  PCU_Max_Doubles(&time_4,1);
460  if (!me) {
461  std::cout<<"\nConstruction time: "<<time_2<<"\n"
462  <<"Problem time: " << time_4<<"\n\n";
463  }
464  //Delete the APF Mesh
465  m->destroyNative();
466  apf::destroyMesh(m);
467  //End communications
468  PCU_Comm_Free();
469 #endif
470 
471  if (me == 0)
472  std::cout << "PASS" << std::endl;
473 
474  return 0;
475 
476 }
477 /*****************************************************************************/
478 /********************************* END MAIN **********************************/
479 /*****************************************************************************/
ColoringProblem sets up coloring problems for the user.
Defines the ColoringProblem class.
#define nParts
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
void solve(bool updateInputData=true)
Direct the problem to create a solution.
MeshAdapter defines the interface for mesh input.
std::bitset< NUM_MODEL_FLAGS > modelFlag_t
int main(int narg, char *arg[])
static ArrayRCP< ArrayRCP< zscalar_t > > weights
Defines the APFMeshAdapter class.
size_t getOwnedList(ArrayView< bool > &isOwner) const
Sets pointer to the ownership of this processes vertices.
size_t getEdgeList(ArrayView< const gno_t > &Ids, ArrayView< input_t > &wgts) const
Sets pointers to this process&#39; hyperedge Ids and their weights.
void PrintGhostMetrics(Zoltan2::HyperGraphModel< Adapter > &mdl)
The StridedData class manages lists of weights or coordinates.
void printTimers() const
Return the communicator passed to the problem.
The user parameters, debug, timing and memory profiling output objects, and error checking methods...
size_t getLocalNumPins() const
Returns the local number of pins.
const PartitioningSolution< Adapter > & getSolution()
Get the solution to the problem.
PartitioningProblem sets up partitioning problems for the user.
Defines the HyperGraphModel interface.
Defines the Environment class.
size_t getPinList(ArrayView< const gno_t > &pinIds, ArrayView< const offset_t > &offsets, ArrayView< input_t > &wgts) const
Sets pointers to this process&#39; pins global Ids based on the centric view given by getCentricView() ...
size_t getLocalNumVertices() const
Returns the number vertices on this process.
Defines the PartitioningProblem class.
A class that computes and returns quality metrics.
HyperGraphModel defines the interface required for hyper graph models.
void solve(bool updateInputData=true)
Direct the problem to create a solution.
size_t getLocalNumOwnedVertices() const
Returns the number vertices on this process that are owned.