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