62 #include "Tpetra_DefaultPlatform.hpp" 65 #include "Teuchos_RCP.hpp" 66 #include "Teuchos_GlobalMPISession.hpp" 67 #include "Teuchos_XMLParameterListHelpers.hpp" 70 #ifdef HAVE_ZOLTAN2_PARMA 81 using Teuchos::ParameterList;
88 typedef Tpetra::DefaultPlatform::DefaultPlatformType
Platform;
90 #ifdef HAVE_ZOLTAN2_PARMA 92 void runTest(RCP<
const Teuchos::Comm<int> >& CommT, apf::Mesh2* m,std::string action,
93 std::string parma_method,
int nParts,
double imbalance, std::string output_title );
100 int main(
int narg,
char *arg[]) {
102 Teuchos::GlobalMPISession mpiSession(&narg, &arg,0);
103 Platform &platform = Tpetra::DefaultPlatform::getDefaultPlatform();
104 RCP<const Teuchos::Comm<int> > CommT = platform.getComm();
106 int me = CommT->getRank();
111 <<
"====================================================================\n" 113 <<
"| Example: Partition APF Mesh |\n" 115 <<
"| Questions? Contact Karen Devine (kddevin@sandia.gov), |\n" 116 <<
"| Erik Boman (egboman@sandia.gov), |\n" 117 <<
"| Siva Rajamanickam (srajama@sandia.gov). |\n" 119 <<
"| Pamgen's website: http://trilinos.sandia.gov/packages/pamgen |\n" 120 <<
"| Zoltan2's website: http://trilinos.sandia.gov/packages/zoltan2 |\n" 121 <<
"| Trilinos website: http://trilinos.sandia.gov |\n" 123 <<
"====================================================================\n";
129 cout <<
"PARALLEL executable \n";
133 cout <<
"SERIAL executable \n";
142 std::string meshFileName(
"4/");
143 std::string modelFileName(
"torus.dmg");
144 std::string action(
"zoltan_hg");
145 std::string parma_method(
"VtxElm");
146 std::string output_loc(
"");
147 int nParts = CommT->getSize();
148 double imbalance=1.1;
151 Teuchos::CommandLineProcessor cmdp (
false,
false);
152 cmdp.setOption(
"meshfile", &meshFileName,
153 "Mesh file with APF specifications (.smb file(s))");
154 cmdp.setOption(
"modelfile", &modelFileName,
155 "Model file with APF specifications (.dmg file)");
156 cmdp.setOption(
"action", &action,
157 "Method to use: mj, scotch, zoltan_rcb, parma or color");
158 cmdp.setOption(
"parma_method", &parma_method,
159 "Method to use: Vertex, Edge, Element, VtxElm, VtxEdgeElm, ElmLtVtx, Ghost, or Shape ");
160 cmdp.setOption(
"nparts", &
nParts,
161 "Number of parts to create");
162 cmdp.setOption(
"imbalance", &imbalance,
163 "Target Imbalance for first partitioner");
164 cmdp.setOption(
"output", &output_loc,
165 "Location of new partitioned apf mesh. Ex: 4/torus.smb");
166 cmdp.parse(narg, arg);
180 #ifdef HAVE_ZOLTAN2_PARMA 182 if (me == 0) cout <<
"Generating mesh ... \n\n";
189 apf::Mesh2* m = apf::loadMdsMesh(modelFileName.c_str(),meshFileName.c_str());
191 runTest(CommT,m,action,parma_method,
nParts,imbalance,
"partition");
193 runTest(CommT,m,
"parma",parma_method,
nParts,imbalance,
"parma");
198 if (output_loc!=
"") {
199 m->writeNative(output_loc.c_str());
203 if (me == 0) cout <<
"Deleting the mesh ... \n\n";
213 std::cout <<
"PASS" << std::endl;
222 #ifdef HAVE_ZOLTAN2_PARMA 224 void runTest(RCP<
const Teuchos::Comm<int> >& CommT, apf::Mesh2* m,std::string action,
225 std::string parma_method,
int nParts,
double imbalance,std::string output_title) {
227 int me = CommT->getRank();
230 std::string primary=
"region";
231 std::string adjacency=
"face";
232 if (m->getDimension()==2) {
236 bool needSecondAdj=
false;
239 if (me == 0) cout <<
"Creating parameter list ... \n\n";
241 Teuchos::ParameterList params(
"test params");
242 params.set(
"timer_output_stream" ,
"std::cout");
244 if (action ==
"mj") {
245 params.set(
"debug_level",
"basic_status");
246 params.set(
"imbalance_tolerance", imbalance);
247 params.set(
"num_global_parts",
nParts);
248 params.set(
"algorithm",
"multijagged");
249 params.set(
"rectilinear",
"yes");
251 else if (action ==
"scotch") {
252 params.set(
"debug_level",
"no_status");
253 params.set(
"imbalance_tolerance", imbalance);
254 params.set(
"num_global_parts",
nParts);
255 params.set(
"partitioning_approach",
"partition");
256 params.set(
"algorithm",
"scotch");
259 else if (action ==
"zoltan_rcb") {
260 params.set(
"debug_level",
"verbose_detailed_status");
261 params.set(
"imbalance_tolerance", imbalance);
262 params.set(
"num_global_parts",
nParts);
263 params.set(
"partitioning_approach",
"partition");
264 params.set(
"algorithm",
"zoltan");
266 else if (action ==
"parma") {
267 params.set(
"debug_level",
"no_status");
268 params.set(
"imbalance_tolerance", imbalance);
269 params.set(
"algorithm",
"parma");
270 Teuchos::ParameterList &pparams = params.sublist(
"parma_parameters",
false);
271 pparams.set(
"parma_method",parma_method);
272 pparams.set(
"step_size",1.1);
273 if (parma_method==
"Ghost") {
274 pparams.set(
"ghost_layers",3);
275 pparams.set(
"ghost_bridge",m->getDimension()-1);
279 else if (action==
"zoltan_hg") {
280 params.set(
"debug_level",
"no_status");
281 params.set(
"imbalance_tolerance", imbalance);
282 params.set(
"algorithm",
"zoltan");
283 params.set(
"num_global_parts",
nParts);
284 Teuchos::ParameterList &zparams = params.sublist(
"zoltan_parameters",
false);
285 zparams.set(
"LB_METHOD",
"HYPERGRAPH");
291 Parma_PrintPtnStats(m,output_title+
"_before");
294 if (me == 0) cout <<
"Creating mesh adapter ... \n\n";
299 double time_1 = PCU_Time();
301 new inputAdapter_t(*CommT, m,primary,adjacency,needSecondAdj);
302 double time_2 = PCU_Time();
306 if (me == 0) cout <<
"Creating partitioning problem ... \n\n";
307 double time_3=PCU_Time();
311 if (me == 0) cout <<
"Calling the partitioner ... \n\n";
317 if (me==0) cout <<
"Applying Solution to Mesh\n\n";
318 apf::Mesh2** new_mesh = &m;
319 ia->applyPartitioningSolution(m,new_mesh,problem.getSolution());
321 double time_4=PCU_Time();
324 RCP<quality_t> metricObject =
325 rcp(
new quality_t(ia, ¶ms, CommT, &problem.getSolution()));
328 metricObject->printMetrics(cout);
332 Parma_PrintPtnStats(m,output_title+
"_after");
337 PCU_Max_Doubles(&time_2,1);
338 PCU_Max_Doubles(&time_4,1);
340 std::cout<<
"\n"<<output_title<<
"Construction time: "<<time_2<<
"\n" 341 <<output_title<<
"Problem time: " << time_4<<
"\n\n";
Defines the ColoringProblem class.
MeshAdapter defines the interface for mesh input.
bool runTest(Adapter &ia, const RCP< const Teuchos::Comm< int > > &comm, std::string hi)
Tpetra::DefaultPlatform::DefaultPlatformType Platform
Defines the APFMeshAdapter class.
int main(int narg, char *arg[])
PartitioningProblem sets up partitioning problems for the user.
Defines the PartitioningProblem class.
A class that computes and returns quality metrics.