63 #include "Tpetra_DefaultPlatform.hpp" 66 #include "Teuchos_RCP.hpp" 67 #include "Teuchos_GlobalMPISession.hpp" 68 #include "Teuchos_XMLParameterListHelpers.hpp" 71 #ifdef HAVE_ZOLTAN2_PARMA 82 using Teuchos::ParameterList;
89 typedef Tpetra::DefaultPlatform::DefaultPlatformType
Platform;
91 #ifdef HAVE_ZOLTAN2_PARMA 93 void runTest(RCP<
const Teuchos::Comm<int> >& CommT, apf::Mesh2* m,std::string action,
94 std::string parma_method,
int nParts,
double imbalance, std::string output_title );
101 int main(
int narg,
char *arg[]) {
103 Teuchos::GlobalMPISession mpiSession(&narg, &arg,0);
104 Platform &platform = Tpetra::DefaultPlatform::getDefaultPlatform();
105 RCP<const Teuchos::Comm<int> > CommT = platform.getComm();
107 int me = CommT->getRank();
112 <<
"====================================================================\n" 114 <<
"| Example: Partition APF Mesh |\n" 116 <<
"| Questions? Contact Karen Devine (kddevin@sandia.gov), |\n" 117 <<
"| Erik Boman (egboman@sandia.gov), |\n" 118 <<
"| Siva Rajamanickam (srajama@sandia.gov). |\n" 120 <<
"| Pamgen's website: http://trilinos.sandia.gov/packages/pamgen |\n" 121 <<
"| Zoltan2's website: http://trilinos.sandia.gov/packages/zoltan2 |\n" 122 <<
"| Trilinos website: http://trilinos.sandia.gov |\n" 124 <<
"====================================================================\n";
130 cout <<
"PARALLEL executable \n";
134 cout <<
"SERIAL executable \n";
143 std::string meshFileName(
"4/");
144 std::string modelFileName(
"torus.dmg");
145 std::string action(
"zoltan_hg");
146 std::string parma_method(
"VtxElm");
147 std::string output_loc(
"");
148 int nParts = CommT->getSize();
149 double imbalance=1.1;
152 Teuchos::CommandLineProcessor cmdp (
false,
false);
153 cmdp.setOption(
"meshfile", &meshFileName,
154 "Mesh file with APF specifications (.smb file(s))");
155 cmdp.setOption(
"modelfile", &modelFileName,
156 "Model file with APF specifications (.dmg file)");
157 cmdp.setOption(
"action", &action,
158 "Method to use: mj, scotch, zoltan_rcb, parma or color");
159 cmdp.setOption(
"parma_method", &parma_method,
160 "Method to use: Vertex, Edge, Element, VtxElm, VtxEdgeElm, ElmLtVtx, Ghost, or Shape ");
161 cmdp.setOption(
"nparts", &nParts,
162 "Number of parts to create");
163 cmdp.setOption(
"imbalance", &imbalance,
164 "Target Imbalance for first partitioner");
165 cmdp.setOption(
"output", &output_loc,
166 "Location of new partitioned apf mesh. Ex: 4/torus.smb");
167 cmdp.parse(narg, arg);
181 #ifdef HAVE_ZOLTAN2_PARMA 183 if (me == 0) cout <<
"Generating mesh ... \n\n";
190 apf::Mesh2* m = apf::loadMdsMesh(modelFileName.c_str(),meshFileName.c_str());
192 runTest(CommT,m,action,parma_method,nParts,imbalance,
"partition");
194 runTest(CommT,m,
"parma",parma_method,nParts,imbalance,
"parma");
199 if (output_loc!=
"") {
200 m->writeNative(output_loc.c_str());
204 if (me == 0) cout <<
"Deleting the mesh ... \n\n";
214 std::cout <<
"PASS" << std::endl;
223 #ifdef HAVE_ZOLTAN2_PARMA 225 void runTest(RCP<
const Teuchos::Comm<int> >& CommT, apf::Mesh2* m,std::string action,
226 std::string parma_method,
int nParts,
double imbalance,std::string output_title) {
228 int me = CommT->getRank();
231 std::string primary=
"region";
232 std::string adjacency=
"face";
233 if (m->getDimension()==2) {
237 bool needSecondAdj=
false;
240 if (me == 0) cout <<
"Creating parameter list ... \n\n";
242 Teuchos::ParameterList params(
"test params");
243 params.set(
"timer_output_stream" ,
"std::cout");
245 if (action ==
"mj") {
246 params.set(
"debug_level",
"basic_status");
247 params.set(
"imbalance_tolerance", imbalance);
248 params.set(
"num_global_parts",
nParts);
249 params.set(
"algorithm",
"multijagged");
250 params.set(
"rectilinear",
"yes");
252 else if (action ==
"scotch") {
253 params.set(
"debug_level",
"no_status");
254 params.set(
"imbalance_tolerance", imbalance);
255 params.set(
"num_global_parts",
nParts);
256 params.set(
"partitioning_approach",
"partition");
257 params.set(
"algorithm",
"scotch");
260 else if (action ==
"zoltan_rcb") {
261 params.set(
"debug_level",
"verbose_detailed_status");
262 params.set(
"imbalance_tolerance", imbalance);
263 params.set(
"num_global_parts",
nParts);
264 params.set(
"partitioning_approach",
"partition");
265 params.set(
"algorithm",
"zoltan");
267 else if (action ==
"parma") {
268 params.set(
"debug_level",
"no_status");
269 params.set(
"imbalance_tolerance", imbalance);
270 params.set(
"algorithm",
"parma");
271 Teuchos::ParameterList &pparams = params.sublist(
"parma_parameters",
false);
272 pparams.set(
"parma_method",parma_method);
273 pparams.set(
"step_size",1.1);
274 if (parma_method==
"Ghost") {
275 pparams.set(
"ghost_layers",3);
276 pparams.set(
"ghost_bridge",m->getDimension()-1);
278 params.set(
"compute_metrics",
"yes");
281 else if (action==
"zoltan_hg") {
282 params.set(
"debug_level",
"no_status");
283 params.set(
"imbalance_tolerance", imbalance);
284 params.set(
"algorithm",
"zoltan");
285 params.set(
"num_global_parts",
nParts);
286 Teuchos::ParameterList &zparams = params.sublist(
"zoltan_parameters",
false);
287 zparams.set(
"LB_METHOD",
"HYPERGRAPH");
293 Parma_PrintPtnStats(m,output_title+
"_before");
296 if (me == 0) cout <<
"Creating mesh adapter ... \n\n";
299 double time_1 = PCU_Time();
300 inputAdapter_t ia(*CommT, m,primary,adjacency,needSecondAdj);
301 double time_2 = PCU_Time();
305 if (me == 0) cout <<
"Creating partitioning problem ... \n\n";
306 double time_3=PCU_Time();
310 if (me == 0) cout <<
"Calling the partitioner ... \n\n";
316 if (me==0) cout <<
"Applying Solution to Mesh\n\n";
317 apf::Mesh2** new_mesh = &m;
318 ia.applyPartitioningSolution(m,new_mesh,problem.
getSolution());
320 double time_4=PCU_Time();
324 Parma_PrintPtnStats(m,output_title+
"_after");
329 PCU_Max_Doubles(&time_2,1);
330 PCU_Max_Doubles(&time_4,1);
332 std::cout<<
"\n"<<output_title<<
"Construction time: "<<time_2<<
"\n" 333 <<output_title<<
"Problem time: " << time_4<<
"\n\n";
Defines the ColoringProblem class.
Tpetra::DefaultPlatform::DefaultPlatformType Platform
Defines the APFMeshAdapter class.
int main(int narg, char *arg[])
void printMetrics(std::ostream &os) const
Print the array of metrics.
const PartitioningSolution< Adapter > & getSolution()
Get the solution to the problem.
PartitioningProblem sets up partitioning problems for the user.
Defines the Environment class.
Defines the PartitioningProblem class.
void solve(bool updateInputData=true)
Direct the problem to create a solution.