Zoltan2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
pamgenMeshAdapterTest.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 
60 
61 // Teuchos includes
62 #include "Teuchos_RCP.hpp"
63 #include "Teuchos_XMLParameterListHelpers.hpp"
64 
65 // Pamgen includes
66 #include "create_inline_mesh.h"
67 
68 using Teuchos::ParameterList;
69 using Teuchos::ArrayRCP;
70 
72 
73 /*****************************************************************************/
74 /******************************** MAIN ***************************************/
75 /*****************************************************************************/
76 
77 int main(int narg, char *arg[]) {
78 
79  Tpetra::ScopeGuard tscope(&narg, &arg);
80  Teuchos::RCP<const Teuchos::Comm<int> > CommT = Tpetra::getDefaultComm();
81 
82  int me = CommT->getRank();
83  int numProcs = CommT->getSize();
84 
85  if (me == 0){
86  std::cout
87  << "====================================================================\n"
88  << "| |\n"
89  << "| Example: Partition Pamgen Hexahedral Mesh |\n"
90  << "| |\n"
91  << "| Questions? Contact Karen Devine (kddevin@sandia.gov), |\n"
92  << "| Erik Boman (egboman@sandia.gov), |\n"
93  << "| Siva Rajamanickam (srajama@sandia.gov). |\n"
94  << "| |\n"
95  << "| Pamgen's website: http://trilinos.sandia.gov/packages/pamgen |\n"
96  << "| Zoltan2's website: http://trilinos.sandia.gov/packages/zoltan2 |\n"
97  << "| Trilinos website: http://trilinos.sandia.gov |\n"
98  << "| |\n"
99  << "====================================================================\n";
100  }
101 
102 
103 #ifdef HAVE_MPI
104  if (me == 0) {
105  std::cout << "PARALLEL executable \n";
106  }
107 #else
108  if (me == 0) {
109  std::cout << "SERIAL executable \n";
110  }
111 #endif
112 
113  /***************************************************************************/
114  /*************************** GET XML INPUTS ********************************/
115  /***************************************************************************/
116 
117  // default values for command-line arguments
118  std::string xmlMeshInFileName("Poisson.xml");
119  std::string action("mj");
120  int nParts = CommT->getSize();
121 
122  // Read run-time options.
123  Teuchos::CommandLineProcessor cmdp (false, false);
124  cmdp.setOption("xmlfile", &xmlMeshInFileName,
125  "XML file with PamGen specifications");
126  cmdp.setOption("action", &action,
127  "Method to use: mj, scotch, zoltan_rcb, zoltan_hg, hg_ghost, "
128  "parma or color");
129  cmdp.setOption("nparts", &nParts,
130  "Number of parts to create");
131  cmdp.parse(narg, arg);
132 
133  // Read xml file into parameter list
134  ParameterList inputMeshList;
135 
136  if(xmlMeshInFileName.length()) {
137  if (me == 0) {
138  std::cout << "\nReading parameter list from the XML file \""
139  <<xmlMeshInFileName<<"\" ...\n\n";
140  }
141  Teuchos::updateParametersFromXmlFile(xmlMeshInFileName,
142  Teuchos::inoutArg(inputMeshList));
143  if (me == 0) {
144  inputMeshList.print(std::cout,2,true,true);
145  std::cout << "\n";
146  }
147  }
148  else {
149  std::cout << "Cannot read input file: " << xmlMeshInFileName << "\n";
150  return 5;
151  }
152 
153  // Get pamgen mesh definition
154  std::string meshInput = Teuchos::getParameter<std::string>(inputMeshList,
155  "meshInput");
156 
157  /***************************************************************************/
158  /********************** GET CELL TOPOLOGY **********************************/
159  /***************************************************************************/
160 
161  // Get dimensions
162  int dim = 3;
163 
164  /***************************************************************************/
165  /***************************** GENERATE MESH *******************************/
166  /***************************************************************************/
167 
168  if (me == 0) std::cout << "Generating mesh ... \n\n";
169 
170  // Generate mesh with Pamgen
171  long long maxInt = 9223372036854775807LL;
172  Create_Pamgen_Mesh(meshInput.c_str(), dim, me, numProcs, maxInt);
173 
174  // Creating mesh adapter
175  if (me == 0) std::cout << "Creating mesh adapter ... \n\n";
176 
177  typedef Zoltan2::PamgenMeshAdapter<basic_user_t> inputAdapter_t;
179 
180  inputAdapter_t *ia = new inputAdapter_t(*CommT, "region");
181  ia->print(me);
182 
183  // Set parameters for partitioning
184  if (me == 0) std::cout << "Creating parameter list ... \n\n";
185 
186  Teuchos::ParameterList params("test params");
187  params.set("timer_output_stream" , "std::cout");
188 
189  bool do_partitioning = false;
190  if (action == "mj") {
191  do_partitioning = true;
192  params.set("debug_level", "basic_status");
193  params.set("imbalance_tolerance", 1.1);
194  params.set("num_global_parts", nParts);
195  params.set("algorithm", "multijagged");
196  params.set("rectilinear", true); // bool parameter
197  }
198  else if (action == "scotch") {
199  do_partitioning = true;
200  params.set("debug_level", "verbose_detailed_status");
201  params.set("imbalance_tolerance", 1.1);
202  params.set("num_global_parts", nParts);
203  params.set("partitioning_approach", "partition");
204  params.set("algorithm", "scotch");
205  }
206  else if (action == "zoltan_rcb") {
207  do_partitioning = true;
208  params.set("debug_level", "verbose_detailed_status");
209  params.set("imbalance_tolerance", 1.1);
210  params.set("num_global_parts", nParts);
211  params.set("partitioning_approach", "partition");
212  params.set("algorithm", "zoltan");
213  }
214  else if (action == "zoltan_hg") {
215  do_partitioning = true;
216  params.set("debug_level", "verbose_detailed_status");
217  params.set("imbalance_tolerance", 1.1);
218  params.set("num_global_parts", nParts);
219  params.set("partitioning_approach", "partition");
220  params.set("algorithm", "zoltan");
221  Teuchos::ParameterList &zparams = params.sublist("zoltan_parameters",false);
222  zparams.set("LB_METHOD","phg");
223  zparams.set("FINAL_OUTPUT", "1");
224  }
225  else if (action=="hg_ghost") {
226  do_partitioning = true;
227  params.set("debug_level", "no_status");
228  params.set("imbalance_tolerance", 1.1);
229  params.set("algorithm", "zoltan");
230  params.set("num_global_parts", nParts);
231  params.set("hypergraph_model_type","ghosting");
232  params.set("ghost_layers",2);
233  Teuchos::ParameterList &zparams = params.sublist("zoltan_parameters",false);
234  zparams.set("LB_METHOD","HYPERGRAPH");
235  zparams.set("LB_APPROACH","PARTITION");
236  zparams.set("PHG_EDGE_SIZE_THRESHOLD", "1.0");
237  }
238 
239  else if (action == "parma") {
240  do_partitioning = true;
241  params.set("debug_level", "basic_status");
242  params.set("imbalance_tolerance", 1.05);
243  params.set("algorithm", "parma");
244  Teuchos::ParameterList &pparams = params.sublist("parma_parameters",false);
245  pparams.set("parma_method","VtxElm");
246  }
247  else if (action=="zoltan_hg") {
248  do_partitioning = true;
249  params.set("debug_level", "no_status");
250  params.set("imbalance_tolerance", 1.1);
251  params.set("algorithm", "zoltan");
252  params.set("num_global_parts", nParts);
253  Teuchos::ParameterList &zparams = params.sublist("zoltan_parameters",false);
254  zparams.set("LB_METHOD","HYPERGRAPH");
255 
256  }
257 
258  else if (action == "color") {
259  params.set("debug_level", "verbose_detailed_status");
260  params.set("debug_output_file", "kdd");
261  params.set("debug_procs", "all");
262  }
263 
264  if(me == 0) std::cout << "Action: " << action << std::endl;
265  // create Partitioning problem
266  if (do_partitioning) {
267  if (me == 0) std::cout << "Creating partitioning problem ... \n\n";
268 
269  Zoltan2::PartitioningProblem<inputAdapter_t> problem(ia, &params, CommT);
270 
271  // call the partitioner
272  if (me == 0) std::cout << "Calling the partitioner ... \n\n";
273 
274  problem.solve();
275 
276  // create metric object
277 
278  Teuchos::RCP<quality_t> metricObject =
279  rcp(new quality_t(ia, &params, CommT, &problem.getSolution()));
280 
281  if (!me) {
282  metricObject->printMetrics(std::cout);
283  }
284  }
285  else {
286  if (me == 0) std::cout << "Creating coloring problem ... \n\n";
287 
288  Zoltan2::ColoringProblem<inputAdapter_t> problem(ia, &params, CommT);
289 
290  // call the partitioner
291  if (me == 0) std::cout << "Calling the coloring algorithm ... \n\n";
292 
293  problem.solve();
294 
295  problem.printTimers();
296  }
297 
298  // delete mesh
299  if (me == 0) std::cout << "Deleting the mesh ... \n\n";
300 
301  Delete_Pamgen_Mesh();
302 
303  if (me == 0)
304  std::cout << "PASS" << std::endl;
305 
306  Kokkos::finalize(); // Test uses directory with kokkos
307 
308  return 0;
309 }
310 /*****************************************************************************/
311 /********************************* END MAIN **********************************/
312 /*****************************************************************************/
ColoringProblem sets up coloring problems for the user.
Defines the ColoringProblem class.
#define nParts
void solve(bool updateInputData=true)
Direct the problem to create a solution.
int main(int narg, char *arg[])
A simple class that can be the User template argument for an InputAdapter.
Defines the PamgenMeshAdapter class.
void printTimers() const
Return the communicator passed to the problem.
const PartitioningSolution< Adapter > & getSolution()
Get the solution to the problem.
PartitioningProblem sets up partitioning problems for the user.
Zoltan2::BasicUserTypes basic_user_t
Defines the PartitioningProblem class.
A class that computes and returns quality metrics.
This class represents a mesh.
void solve(bool updateInputData=true)
Direct the problem to create a solution.