63 #include "Tpetra_DefaultPlatform.hpp" 66 #include "Teuchos_RCP.hpp" 67 #include "Teuchos_GlobalMPISession.hpp" 68 #include "Teuchos_XMLParameterListHelpers.hpp" 71 #include "create_inline_mesh.h" 74 using Teuchos::ParameterList;
81 typedef Tpetra::DefaultPlatform::DefaultPlatformType
Platform;
90 int main(
int narg,
char *arg[]) {
92 Teuchos::GlobalMPISession mpiSession(&narg, &arg,0);
93 Platform &platform = Tpetra::DefaultPlatform::getDefaultPlatform();
94 RCP<const Teuchos::Comm<int> > CommT = platform.getComm();
96 int me = CommT->getRank();
97 int numProcs = CommT->getSize();
101 <<
"====================================================================\n" 103 <<
"| Example: Partition Pamgen Hexahedral Mesh |\n" 105 <<
"| Questions? Contact Karen Devine (kddevin@sandia.gov), |\n" 106 <<
"| Erik Boman (egboman@sandia.gov), |\n" 107 <<
"| Siva Rajamanickam (srajama@sandia.gov). |\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" 113 <<
"====================================================================\n";
119 cout <<
"PARALLEL executable \n";
123 cout <<
"SERIAL executable \n";
132 std::string xmlMeshInFileName(
"Poisson.xml");
133 std::string action(
"mj");
134 int nParts = CommT->getSize();
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, " 143 cmdp.setOption(
"nparts", &nParts,
144 "Number of parts to create");
145 cmdp.parse(narg, arg);
148 ParameterList inputMeshList;
150 if(xmlMeshInFileName.length()) {
152 cout <<
"\nReading parameter list from the XML file \"" 153 <<xmlMeshInFileName<<
"\" ...\n\n";
155 Teuchos::updateParametersFromXmlFile(xmlMeshInFileName,
156 Teuchos::inoutArg(inputMeshList));
158 inputMeshList.print(cout,2,
true,
true);
163 cout <<
"Cannot read input file: " << xmlMeshInFileName <<
"\n";
168 std::string meshInput = Teuchos::getParameter<std::string>(inputMeshList,
182 if (me == 0) cout <<
"Generating mesh ... \n\n";
185 long long maxInt = 9223372036854775807LL;
186 Create_Pamgen_Mesh(meshInput.c_str(), dim, me, numProcs, maxInt);
189 if (me == 0) cout <<
"Creating mesh adapter ... \n\n";
193 inputAdapter_t ia(*CommT,
"region");
197 if (me == 0) cout <<
"Creating parameter list ... \n\n";
199 Teuchos::ParameterList params(
"test params");
200 params.set(
"timer_output_stream" ,
"std::cout");
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");
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");
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");
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");
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");
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");
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");
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");
279 if(me == 0) cout <<
"Action: " << action << endl;
281 if (do_partitioning) {
282 if (me == 0) cout <<
"Creating partitioning problem ... \n\n";
287 if (me == 0) cout <<
"Calling the partitioner ... \n\n";
294 if (action ==
"scotch")
299 if (me == 0) cout <<
"Creating coloring problem ... \n\n";
304 if (me == 0) cout <<
"Calling the coloring algorithm ... \n\n";
312 if (me == 0) cout <<
"Deleting the mesh ... \n\n";
314 Delete_Pamgen_Mesh();
317 std::cout <<
"PASS" << std::endl;
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
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.
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