64 #include "Teuchos_RCP.hpp"
65 #include "Teuchos_XMLParameterListHelpers.hpp"
66 #include "Teuchos_Hashtable.hpp"
69 #ifdef HAVE_ZOLTAN2_PARMA
81 using Teuchos::ParameterList;
83 using Teuchos::ArrayView;
86 template <
typename Adapter>
90 typedef typename Adapter::scalar_t scalar_t;
91 typedef typename Adapter::offset_t offset_t;
94 ArrayView<const gno_t> Ids;
95 ArrayView<input_t> wgts;
97 ArrayView<bool> isOwner;
100 ArrayView<const gno_t> pins;
101 ArrayView<const offset_t> offsets;
104 std::set<gno_t> gids;
109 std::set<gno_t> ghosts;
113 if (gids.find(pin)==gids.end()) {
115 if (ghosts.find(pin)==ghosts.end())
119 std::cout<<
"[METRIC] " << PCU_Comm_Self() <<
" Total number of ghosts in the hypergraph: " << num_ghosts <<
"\n"
120 <<
"[METRIC] " << PCU_Comm_Self() <<
" Number of unique ghosts: " << ghosts.size() <<
"\n";
121 gno_t unique_ghosts =ghosts.size();
122 gno_t owned_and_ghosts =unique_ghosts+numOwned;
123 gno_t max_o_and_g,min_o_and_g;
124 gno_t max_ghosts,max_u_ghosts;
125 gno_t min_ghosts,min_u_ghosts;
126 max_ghosts = min_ghosts = num_ghosts;
127 max_u_ghosts = min_u_ghosts = unique_ghosts;
128 max_o_and_g = min_o_and_g = owned_and_ghosts;
129 double avg_ghosts,avg_u_ghosts,avg_o_and_g;
130 PCU_Add_Ints(&num_ghosts,1);
131 PCU_Add_Ints(&unique_ghosts,1);
132 PCU_Add_Ints(&owned_and_ghosts,1);
133 PCU_Max_Ints(&max_ghosts,1);
134 PCU_Max_Ints(&max_u_ghosts,1);
135 PCU_Max_Ints(&max_o_and_g,1);
136 PCU_Min_Ints(&min_ghosts,1);
137 PCU_Min_Ints(&min_u_ghosts,1);
138 PCU_Min_Ints(&min_o_and_g,1);
139 avg_ghosts = num_ghosts*1.0/PCU_Comm_Peers();
140 avg_u_ghosts = unique_ghosts*1.0/PCU_Comm_Peers();
141 avg_o_and_g = owned_and_ghosts*1.0/PCU_Comm_Peers();
142 if (!PCU_Comm_Self())
143 std::cout<<
"[METRIC] Global ghosts in the hypergraph (tot max min avg imb): "
144 << num_ghosts<<
" "<<max_ghosts<<
" "<<min_ghosts<<
" "<<avg_ghosts<<
" "
145 <<max_ghosts/avg_ghosts <<
"\n"
146 <<
"[METRIC] Global unique ghosts (tot max min avg imb): "
147 << unique_ghosts<<
" "<<max_u_ghosts<<
" "<<min_u_ghosts<<
" "<<avg_u_ghosts<<
" "
148 <<max_u_ghosts/avg_u_ghosts <<
"\n"
149 <<
"[METRIC] Global owned and ghosts (tot max min avg imb): "
150 << owned_and_ghosts<<
" "<<max_o_and_g<<
" "<<min_o_and_g<<
" "<<avg_o_and_g<<
" "
151 <<max_o_and_g/avg_o_and_g <<
"\n";
159 int main(
int narg,
char *arg[]) {
161 Tpetra::ScopeGuard tscope(&narg, &arg);
162 Teuchos::RCP<const Teuchos::Comm<int> > CommT = Tpetra::getDefaultComm();
164 int me = CommT->getRank();
169 <<
"====================================================================\n"
171 <<
"| Example: Partition APF Mesh |\n"
173 <<
"| Questions? Contact Karen Devine (kddevin@sandia.gov), |\n"
174 <<
"| Erik Boman (egboman@sandia.gov), |\n"
175 <<
"| Siva Rajamanickam (srajama@sandia.gov). |\n"
177 <<
"| Zoltan2's website: http://trilinos.sandia.gov/packages/zoltan2 |\n"
178 <<
"| Trilinos website: http://trilinos.sandia.gov |\n"
180 <<
"====================================================================\n";
186 std::cout <<
"PARALLEL executable \n";
190 std::cout <<
"SERIAL executable \n";
199 std::string meshFileName(
"4/");
200 std::string modelFileName(
"torus.dmg");
201 std::string action(
"parma");
202 std::string parma_method(
"VtxElm");
203 std::string output_loc(
"");
204 int nParts = CommT->getSize();
205 double imbalance = 1.1;
207 int ghost_metric=
false;
209 Teuchos::CommandLineProcessor cmdp (
false,
false);
210 cmdp.setOption(
"meshfile", &meshFileName,
211 "Mesh file with APF specifications (.smb file(s))");
212 cmdp.setOption(
"modelfile", &modelFileName,
213 "Model file with APF specifications (.dmg file)");
214 cmdp.setOption(
"action", &action,
215 "Method to use: mj, scotch, zoltan_rcb, parma or color");
216 cmdp.setOption(
"parma_method", &parma_method,
217 "Method to use: Vertex, Element, VtxElm, VtxEdgeElm, Ghost, Shape, or Centroid ");
218 cmdp.setOption(
"nparts", &nParts,
219 "Number of parts to create");
220 cmdp.setOption(
"imbalance", &imbalance,
221 "Target imbalance for the partitioning method");
222 cmdp.setOption(
"output", &output_loc,
223 "Location of new partitioned apf mesh. Ex: 4/torus.smb");
224 cmdp.setOption(
"layers", &layers,
225 "Number of layers for ghosting");
226 cmdp.setOption(
"ghost_metric", &ghost_metric,
227 "0 does not compute ghost metric otherwise compute both before and after");
228 cmdp.parse(narg, arg);
241 #ifdef HAVE_ZOLTAN2_PARMA
243 if (me == 0) std::cout <<
"Generating mesh ... \n\n";
250 apf::Mesh2* m = apf::loadMdsMesh(modelFileName.c_str(),meshFileName.c_str());
254 std::string primary=
"region";
255 std::string adjacency=
"face";
256 if (m->getDimension()==2) {
260 bool needSecondAdj=
false;
263 if (me == 0) std::cout <<
"Creating parameter list ... \n\n";
265 Teuchos::ParameterList params(
"test params");
266 params.set(
"timer_output_stream" ,
"std::cout");
268 bool do_partitioning =
false;
269 if (action ==
"mj") {
270 do_partitioning =
true;
271 params.set(
"debug_level",
"basic_status");
272 params.set(
"imbalance_tolerance", imbalance);
273 params.set(
"num_global_parts", nParts);
274 params.set(
"algorithm",
"multijagged");
275 params.set(
"rectilinear",
"yes");
277 else if (action ==
"scotch") {
278 do_partitioning =
true;
279 params.set(
"debug_level",
"no_status");
280 params.set(
"imbalance_tolerance", imbalance);
281 params.set(
"num_global_parts", nParts);
282 params.set(
"partitioning_approach",
"partition");
283 params.set(
"objects_to_partition",
"mesh_elements");
284 params.set(
"algorithm",
"scotch");
287 else if (action ==
"zoltan_rcb") {
288 do_partitioning =
true;
289 params.set(
"debug_level",
"verbose_detailed_status");
290 params.set(
"imbalance_tolerance", imbalance);
291 params.set(
"num_global_parts", nParts);
292 params.set(
"partitioning_approach",
"partition");
293 params.set(
"algorithm",
"zoltan");
295 else if (action ==
"parma") {
296 do_partitioning =
true;
297 params.set(
"debug_level",
"no_status");
298 params.set(
"imbalance_tolerance", imbalance);
299 params.set(
"algorithm",
"parma");
300 Teuchos::ParameterList &pparams = params.sublist(
"parma_parameters",
false);
301 pparams.set(
"parma_method",parma_method);
302 pparams.set(
"step_size",1.1);
303 if (parma_method==
"Ghost") {
304 pparams.set(
"ghost_layers",layers);
305 pparams.set(
"ghost_bridge",m->getDimension()-1);
309 else if (action==
"zoltan_hg") {
310 do_partitioning =
true;
311 params.set(
"debug_level",
"no_status");
312 params.set(
"imbalance_tolerance", imbalance);
313 params.set(
"algorithm",
"zoltan");
314 params.set(
"num_global_parts", nParts);
315 Teuchos::ParameterList &zparams = params.sublist(
"zoltan_parameters",
false);
316 zparams.set(
"LB_METHOD",
"HYPERGRAPH");
317 zparams.set(
"LB_APPROACH",
"PARTITION");
321 else if (action==
"hg_ghost") {
322 do_partitioning =
true;
323 params.set(
"debug_level",
"no_status");
324 params.set(
"imbalance_tolerance", imbalance);
325 params.set(
"algorithm",
"zoltan");
326 params.set(
"num_global_parts", nParts);
327 params.set(
"hypergraph_model_type",
"ghosting");
328 params.set(
"ghost_layers",layers);
329 Teuchos::ParameterList &zparams = params.sublist(
"zoltan_parameters",
false);
330 zparams.set(
"LB_METHOD",
"HYPERGRAPH");
331 zparams.set(
"LB_APPROACH",
"PARTITION");
332 zparams.set(
"PHG_EDGE_SIZE_THRESHOLD",
"1.0");
337 else if (action ==
"color") {
338 params.set(
"debug_level",
"verbose_detailed_status");
339 params.set(
"debug_output_file",
"kdd");
340 params.set(
"debug_procs",
"all");
342 Parma_PrintPtnStats(m,
"before");
345 if (me == 0) std::cout <<
"Creating mesh adapter ... \n\n";
350 double time_1=PCU_Time();
352 new inputAdapter_t(*CommT, m,primary,adjacency,needSecondAdj);
353 double time_2=PCU_Time();
356 inputAdapter_t::scalar_t* arr =
357 new inputAdapter_t::scalar_t[ia->getLocalNumOf(ia->getPrimaryEntityType())];
358 for (
size_t i=0;i<ia->getLocalNumOf(ia->getPrimaryEntityType());i++) {
359 arr[i]=PCU_Comm_Self()+1;
362 const inputAdapter_t::scalar_t*
weights=arr;
363 ia->setWeights(ia->getPrimaryEntityType(),
weights,1);
367 const baseMeshAdapter_t *base_ia =
dynamic_cast<const baseMeshAdapter_t*
>(ia);
369 RCP<Zoltan2::Environment> env;
375 RCP<const Zoltan2::Environment> envConst = Teuchos::rcp_const_cast<
const Zoltan2::Environment>(env);
377 RCP<const baseMeshAdapter_t> baseInputAdapter_(base_ia,
false);
384 double time_3 = PCU_Time();
385 if (do_partitioning) {
386 if (me == 0) std::cout <<
"Creating partitioning problem ... \n\n";
391 if (me == 0) std::cout <<
"Calling the partitioner ... \n\n";
397 if (me==0) std::cout <<
"Applying Solution to Mesh\n\n";
398 apf::Mesh2** new_mesh = &m;
399 ia->applyPartitioningSolution(m,new_mesh,problem.
getSolution());
402 RCP<quality_t> metricObject =
406 metricObject->printMetrics(std::cout);
410 if (me == 0) std::cout <<
"Creating coloring problem ... \n\n";
415 if (me == 0) std::cout <<
"Calling the coloring algorithm ... \n\n";
424 double time_4=PCU_Time();
432 inputAdapter_t ia2(*CommT, m,primary,adjacency,
true);
433 const baseMeshAdapter_t *base_ia =
dynamic_cast<const baseMeshAdapter_t*
>(&ia2);
436 RCP<Zoltan2::Environment> env;
441 RCP<const Zoltan2::Environment> envConst = Teuchos::rcp_const_cast<
const Zoltan2::Environment>(env);
442 RCP<const baseMeshAdapter_t> baseInputAdapter_(base_ia,
false);
450 if (output_loc!=
"") {
451 m->writeNative(output_loc.c_str());
455 if (me == 0) std::cout <<
"Deleting the mesh ... \n\n";
458 PCU_Max_Doubles(&time_2,1);
459 PCU_Max_Doubles(&time_4,1);
461 std::cout<<
"\nConstruction time: "<<time_2<<
"\n"
462 <<
"Problem time: " << time_4<<
"\n\n";
472 std::cout <<
"PASS" << std::endl;
ColoringProblem sets up coloring problems for the user.
Defines the ColoringProblem class.
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
void solve(bool updateInputData=true)
Direct the problem to create a solution.
MeshAdapter defines the interface for mesh input.
std::bitset< NUM_MODEL_FLAGS > modelFlag_t
map_t::global_ordinal_type gno_t
Zoltan2::EvaluatePartition< matrixAdapter_t > quality_t
int main(int narg, char **arg)
Defines the APFMeshAdapter class.
size_t getOwnedList(ArrayView< bool > &isOwner) const
Sets pointer to the ownership of this processes vertices.
size_t getEdgeList(ArrayView< const gno_t > &Ids, ArrayView< input_t > &wgts) const
Sets pointers to this process' hyperedge Ids and their weights.
void PrintGhostMetrics(Zoltan2::HyperGraphModel< Adapter > &mdl)
The StridedData class manages lists of weights or coordinates.
map_t::local_ordinal_type lno_t
void printTimers() const
Return the communicator passed to the problem.
The user parameters, debug, timing and memory profiling output objects, and error checking methods...
size_t getLocalNumPins() const
Returns the local number of pins.
const PartitioningSolution< Adapter > & getSolution()
Get the solution to the problem.
PartitioningProblem sets up partitioning problems for the user.
Defines the HyperGraphModel interface.
Defines the Environment class.
size_t getPinList(ArrayView< const gno_t > &pinIds, ArrayView< const offset_t > &offsets, ArrayView< input_t > &wgts) const
Sets pointers to this process' pins global Ids based on the centric view given by getCentricView() ...
size_t getLocalNumVertices() const
Returns the number vertices on this process.
Defines the PartitioningProblem class.
A class that computes and returns quality metrics.
HyperGraphModel defines the interface required for hyper graph models.
void solve(bool updateInputData=true)
Direct the problem to create a solution.
size_t getLocalNumOwnedVertices() const
Returns the number vertices on this process that are owned.