62 #include "Tpetra_DefaultPlatform.hpp" 65 #include "Teuchos_RCP.hpp" 66 #include "Teuchos_GlobalMPISession.hpp" 67 #include "Teuchos_XMLParameterListHelpers.hpp" 70 #include "create_inline_mesh.h" 73 using Teuchos::ParameterList;
75 using Teuchos::ArrayRCP;
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 typedef inputAdapter_t::part_t part_t;
196 inputAdapter_t *ia =
new inputAdapter_t(*CommT,
"region");
200 if (me == 0) cout <<
"Creating parameter list ... \n\n";
202 Teuchos::ParameterList params(
"test params");
203 params.set(
"timer_output_stream" ,
"std::cout");
205 bool do_partitioning =
false;
206 if (action ==
"mj") {
207 do_partitioning =
true;
208 params.set(
"debug_level",
"basic_status");
209 params.set(
"imbalance_tolerance", 1.1);
210 params.set(
"num_global_parts",
nParts);
211 params.set(
"algorithm",
"multijagged");
212 params.set(
"rectilinear",
true);
214 else if (action ==
"scotch") {
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",
"scotch");
222 else if (action ==
"zoltan_rcb") {
223 do_partitioning =
true;
224 params.set(
"debug_level",
"verbose_detailed_status");
225 params.set(
"imbalance_tolerance", 1.1);
226 params.set(
"num_global_parts",
nParts);
227 params.set(
"partitioning_approach",
"partition");
228 params.set(
"algorithm",
"zoltan");
230 else if (action ==
"zoltan_hg") {
231 do_partitioning =
true;
232 params.set(
"debug_level",
"verbose_detailed_status");
233 params.set(
"imbalance_tolerance", 1.1);
234 params.set(
"num_global_parts",
nParts);
235 params.set(
"partitioning_approach",
"partition");
236 params.set(
"algorithm",
"zoltan");
237 Teuchos::ParameterList &zparams = params.sublist(
"zoltan_parameters",
false);
238 zparams.set(
"LB_METHOD",
"phg");
239 zparams.set(
"FINAL_OUTPUT",
"1");
241 else if (action==
"hg_ghost") {
242 do_partitioning =
true;
243 params.set(
"debug_level",
"no_status");
244 params.set(
"imbalance_tolerance", 1.1);
245 params.set(
"algorithm",
"zoltan");
246 params.set(
"num_global_parts",
nParts);
247 params.set(
"hypergraph_model_type",
"ghosting");
248 params.set(
"ghost_layers",2);
249 Teuchos::ParameterList &zparams = params.sublist(
"zoltan_parameters",
false);
250 zparams.set(
"LB_METHOD",
"HYPERGRAPH");
251 zparams.set(
"LB_APPROACH",
"PARTITION");
252 zparams.set(
"PHG_EDGE_SIZE_THRESHOLD",
"1.0");
255 else if (action ==
"parma") {
256 do_partitioning =
true;
257 params.set(
"debug_level",
"basic_status");
258 params.set(
"imbalance_tolerance", 1.05);
259 params.set(
"algorithm",
"parma");
260 Teuchos::ParameterList &pparams = params.sublist(
"parma_parameters",
false);
261 pparams.set(
"parma_method",
"VtxElm");
263 else if (action==
"zoltan_hg") {
264 do_partitioning =
true;
265 params.set(
"debug_level",
"no_status");
266 params.set(
"imbalance_tolerance", 1.1);
267 params.set(
"algorithm",
"zoltan");
268 params.set(
"num_global_parts",
nParts);
269 Teuchos::ParameterList &zparams = params.sublist(
"zoltan_parameters",
false);
270 zparams.set(
"LB_METHOD",
"HYPERGRAPH");
274 else if (action ==
"color") {
275 params.set(
"debug_level",
"verbose_detailed_status");
276 params.set(
"debug_output_file",
"kdd");
277 params.set(
"debug_procs",
"all");
280 if(me == 0) cout <<
"Action: " << action << endl;
282 if (do_partitioning) {
283 if (me == 0) cout <<
"Creating partitioning problem ... \n\n";
288 if (me == 0) cout <<
"Calling the partitioner ... \n\n";
294 RCP<quality_t> metricObject =
295 rcp(
new quality_t(ia, ¶ms, CommT, &problem.
getSolution()));
298 metricObject->printMetrics(cout);
302 if (me == 0) cout <<
"Creating coloring problem ... \n\n";
307 if (me == 0) cout <<
"Calling the coloring algorithm ... \n\n";
315 if (me == 0) cout <<
"Deleting the mesh ... \n\n";
317 Delete_Pamgen_Mesh();
320 std::cout <<
"PASS" << std::endl;
Zoltan2::BaseAdapter< userTypes_t > base_adapter_t
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 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.
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.
Tpetra::MultiVector< double > tMVector_t