Zoltan2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
TaskMappingTest.cpp
Go to the documentation of this file.
1 
5 
6 #include <string>
7 
8 #include <Teuchos_RCP.hpp>
9 #include <Teuchos_Array.hpp>
10 #include <Teuchos_ParameterList.hpp>
11 #include "Teuchos_XMLParameterListHelpers.hpp"
12 
13 #include "Tpetra_MultiVector.hpp"
14 #include <Tpetra_CrsGraph.hpp>
15 #include <Tpetra_Map.hpp>
16 
19 
20 typedef Tpetra::CrsGraph<zlno_t, zgno_t, znode_t> tcrsGraph_t;
21 typedef Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t> tMVector_t;
23 
24 
25 int main(int narg, char *arg[]) {
26 
27  Tpetra::ScopeGuard tscope(&narg, &arg);
28  Teuchos::RCP<const Teuchos::Comm<int> > tcomm = Tpetra::getDefaultComm();
29 
31 
32  int nx = 2, ny = 2, nz = 2;
33  for (int i = 1 ; i < narg ; ++i) {
34  if (0 == strcasecmp( arg[i] , "NX")) {
35  nx = atoi( arg[++i] );
36  }
37  else if (0 == strcasecmp( arg[i] , "NY")) {
38  ny = atoi( arg[++i] );
39  }
40  else if (0 == strcasecmp( arg[i] , "NZ")) {
41  nz = atoi( arg[++i] );
42  }
43  else{
44  std::cerr << "Unrecognized command line argument #"
45  << i << ": " << arg[i] << std::endl ;
46  return 1;
47  }
48  }
49 
50  try{
51 
52  int rank = tcomm->getRank();
53  part_t numProcs = tcomm->getSize();
54 
55  int coordDim = 3;
56  zgno_t numGlobalTasks = nx*ny*nz;
57 
58  zgno_t myTasks = numGlobalTasks / numProcs;
59  zgno_t taskLeftOver = numGlobalTasks % numProcs;
60  if (rank < taskLeftOver ) ++myTasks;
61 
62  zgno_t myTaskBegin = (numGlobalTasks / numProcs) * rank;
63  myTaskBegin += (taskLeftOver < rank ? taskLeftOver : rank);
64  zgno_t myTaskEnd = myTaskBegin + myTasks;
65 
66  zscalar_t **partCenters = NULL;
67  partCenters = new zscalar_t * [coordDim];
68  for(int i = 0; i < coordDim; ++i) {
69  partCenters[i] = new zscalar_t[myTasks];
70  }
71 
72  zgno_t *task_gnos = new zgno_t [myTasks];
73  zlno_t *task_communication_xadj_ = new zlno_t [myTasks + 1];
74  zgno_t *task_communication_adj_ = new zgno_t [myTasks * 6];
75 
76  zlno_t prevNCount = 0;
77  task_communication_xadj_[0] = 0;
78  for (zlno_t i = myTaskBegin; i < myTaskEnd; ++i) {
79  task_gnos[i - myTaskBegin] = i;
80 
81  zlno_t x = i % nx;
82  zlno_t y = (i / (nx)) % ny;
83  zlno_t z = (i / (nx)) / ny;
84  partCenters[0][i - myTaskBegin] = x;
85  partCenters[1][i - myTaskBegin] = y;
86  partCenters[2][i - myTaskBegin] = z;
87 
88  if (x > 0) {
89  task_communication_adj_[prevNCount++] = i - 1;
90  }
91  if (x < nx - 1) {
92  task_communication_adj_[prevNCount++] = i + 1;
93  }
94  if (y > 0) {
95  task_communication_adj_[prevNCount++] = i - nx;
96  }
97  if (y < ny - 1) {
98  task_communication_adj_[prevNCount++] = i + nx;
99  }
100  if (z > 0) {
101  task_communication_adj_[prevNCount++] = i - nx * ny;
102  }
103  if (z < nz - 1) {
104  task_communication_adj_[prevNCount++] = i + nx * ny;
105  }
106  task_communication_xadj_[i + 1 - myTaskBegin] = prevNCount;
107  }
108  using namespace Teuchos;
109  RCP<my_adapter_t> ia;
110  typedef Tpetra::Map<>::node_type mytest_znode_t;
111  typedef Tpetra::Map<zlno_t, zgno_t, mytest_znode_t> map_t;
112  RCP<const map_t> map =
113  rcp(new map_t (numGlobalTasks, myTasks, 0, tcomm));
114 
115  Teuchos::Array<size_t> adjPerTask(myTasks);
116  for (zlno_t lclRow = 0; lclRow < myTasks; lclRow++)
117  adjPerTask[lclRow] = task_communication_xadj_[lclRow+1]
118  - task_communication_xadj_[lclRow];
119  RCP<tcrsGraph_t> TpetraCrsGraph(new tcrsGraph_t (map, adjPerTask()));
120 
121  for (zlno_t lclRow = 0; lclRow < myTasks; ++lclRow) {
122  const zgno_t gblRow = map->getGlobalElement (lclRow);
123  zgno_t begin = task_communication_xadj_[lclRow];
124  zgno_t end = task_communication_xadj_[lclRow + 1];
125  const ArrayView< const zgno_t > indices(task_communication_adj_ + begin,
126  end - begin);
127  TpetraCrsGraph->insertGlobalIndices(gblRow, indices);
128  }
129  TpetraCrsGraph->fillComplete ();
130  RCP<const tcrsGraph_t> const_data =
131  rcp_const_cast<const tcrsGraph_t>(TpetraCrsGraph);
132 
133  ia = RCP<my_adapter_t> (new my_adapter_t(const_data));
134 
135  const int coord_dim = 3;
136  Teuchos::Array<Teuchos::ArrayView<const zscalar_t> > coordView(coord_dim);
137 
138  if(myTasks > 0) {
139  Teuchos::ArrayView<const zscalar_t> a(partCenters[0], myTasks);
140  coordView[0] = a;
141  Teuchos::ArrayView<const zscalar_t> b(partCenters[1], myTasks);
142  coordView[1] = b;
143  Teuchos::ArrayView<const zscalar_t> c(partCenters[2], myTasks);
144  coordView[2] = c;
145  }
146  else {
147  Teuchos::ArrayView<const zscalar_t> a;
148  coordView[0] = a;
149  coordView[1] = a;
150  coordView[2] = a;
151  }
152  RCP<tMVector_t> coords(new tMVector_t(map, coordView.view(0, coord_dim),
153  coord_dim));//= set multivector;
154  RCP<const tMVector_t> const_coords =
155  rcp_const_cast<const tMVector_t>(coords);
156  RCP <Zoltan2::XpetraMultiVectorAdapter<tMVector_t> > adapter(
158  ia->setCoordinateInput(adapter.getRawPtr());
159 
160 // return ia;
161 
162  // Create input adapter
163 // RCP<my_adapter_t> ia = create_problem(tcomm, nx, ny, nz);
164 
165  // Create partitioning problem
166  // xpetra_graph problem type
167  typedef Zoltan2::PartitioningProblem<my_adapter_t> xcrsGraph_problem_t;
169  ParameterList zoltan2_parameters;
170  zoltan2_parameters.set("compute_metrics", true); // bool parameter
171  zoltan2_parameters.set("imbalance_tolerance", 1.0);
172  zoltan2_parameters.set("num_global_parts", tcomm->getSize());
173  zoltan2_parameters.set("algorithm", "multijagged");
174  zoltan2_parameters.set("mj_keep_part_boxes", false); // bool parameter
175  zoltan2_parameters.set("mj_recursion_depth", 3);
176 
177  RCP<xcrsGraph_problem_t> partition_problem;
178  partition_problem =
179  RCP<xcrsGraph_problem_t> (new xcrsGraph_problem_t(
180  ia.getRawPtr(),&zoltan2_parameters,tcomm));
181 
182  // Solve the partitioning problem.
183  partition_problem->solve();
184  tcomm->barrier();
185  RCP<const Zoltan2::Environment> env = partition_problem->getEnvironment();
186 
187  RCP<quality_t>metricObject =
188  rcp(new quality_t(ia.getRawPtr(), &zoltan2_parameters, tcomm,
189  &partition_problem->getSolution()));
190 
191  if (tcomm->getRank() == 0) {
192  metricObject->printMetrics(std::cout);
193  }
194  partition_problem->printTimers();
195 
196  part_t *proc_to_task_xadj_ = NULL;
197  part_t *proc_to_task_adj_ = NULL;
198 
199  // Create the zoltan2 machine representation object
201 
202  // Create the mapper and map the partitioning solution.
204  tcomm,
205  Teuchos::rcpFromRef(mach),
206  ia,
207  rcpFromRef(partition_problem->getSolution()),
208  env);
209 
210  // Get the results and print
211  ctm.getProcTask(proc_to_task_xadj_, proc_to_task_adj_);
212 // part_t numProcs = tcomm->getSize();
213  if (tcomm->getRank() == 0) {
214  for (part_t i = 0; i < numProcs; ++i) {
215  std::cout << "\nProc i:" << i << " ";
216  for (part_t j = proc_to_task_xadj_[i];
217  j < proc_to_task_xadj_[i + 1]; ++j) {
218  std::cout << " " << proc_to_task_adj_[j];
219  }
220  }
221  std::cout << std::endl;
222  }
223 
224  // Below is to calculate the result hops. this uses the global graph
225  // also this is used for debug, as the hops are also calculated in mapper.
226  {
227  zlno_t prevNCount_tmp = 0;
228  zgno_t *task_communication_adj_tmp = new zgno_t [numGlobalTasks * 6];
229  zlno_t *task_communication_xadj_tmp = new zlno_t [numGlobalTasks + 1];
230  task_communication_xadj_tmp[0] = 0;
231 
232  for (zlno_t i = 0; i < numGlobalTasks; ++i) {
233  zlno_t x = i % nx;
234  zlno_t y = (i / (nx)) % ny;
235  zlno_t z = (i / (nx)) / ny;
236 
237  if (x > 0) {
238  task_communication_adj_tmp[prevNCount_tmp++] = i - 1;
239  }
240  if (x < nx - 1) {
241  task_communication_adj_tmp[prevNCount_tmp++] = i + 1;
242  }
243  if (y > 0) {
244  task_communication_adj_tmp[prevNCount_tmp++] = i - nx;
245  }
246  if (y < ny - 1) {
247  task_communication_adj_tmp[prevNCount_tmp++] = i + nx;
248  }
249  if (z > 0) {
250  task_communication_adj_tmp[prevNCount_tmp++] = i - nx * ny;
251  }
252  if (z < nz - 1) {
253  task_communication_adj_tmp[prevNCount_tmp++] = i + nx * ny;
254  }
255  task_communication_xadj_tmp[i + 1] = prevNCount_tmp;
256  }
257 
258  int mach_coord_dim = mach.getMachineDim();
259  std::vector <int> mach_extent(mach_coord_dim);
260  mach.getMachineExtent(&(mach_extent[0]));
261 
262  std::vector <part_t> all_parts(numGlobalTasks), copy(numGlobalTasks, 0);
263 
264  const part_t *parts =
265  partition_problem->getSolution().getPartListView();
266 
267 // typedef Tpetra::Map<>::node_type mytest_znode_t;
268 // typedef Tpetra::Map<zlno_t, zgno_t, mytest_znode_t> map_t;
269 // RCP<const map_t> map =
270 // rcp(new map_t (numGlobalTasks, myTasks, 0, tcomm));
271 
272  for (part_t i = 0; i < myTasks; ++i) {
273  zgno_t g = map->getGlobalElement(i);
274  copy[g] = parts[i];
275  }
276 
277  reduceAll<int, part_t>(
278  *tcomm,
279  Teuchos::REDUCE_SUM,
280  numGlobalTasks,
281  &(copy[0]),
282  &(all_parts[0])
283  );
284 
285  zscalar_t **proc_coords;
286  mach.getAllMachineCoordinatesView(proc_coords);
287  part_t hops=0;
288  part_t hops2 = 0;
289  int *machine_extent = new int [mach_coord_dim];
290  bool *machine_extent_wrap_around = new bool[mach_coord_dim];
291 
292  mach.getMachineExtent(machine_extent);
293  mach.getMachineExtentWrapArounds(machine_extent_wrap_around);
294 
295  for (zlno_t i = 0; i < numGlobalTasks; ++i) {
296  zlno_t b = task_communication_xadj_tmp[i];
297  zlno_t e = task_communication_xadj_tmp[i + 1];
298 
299  part_t procId1 = ctm.getAssignedProcForTask(all_parts[i]);
300 
301  for (zlno_t j = b; j < e; ++j) {
302  zgno_t n = task_communication_adj_tmp[j];
303  part_t procId2 = ctm.getAssignedProcForTask(all_parts[n]);
304 
305  zscalar_t distance2 = 0;
306  mach.getHopCount(procId1, procId2, distance2);
307 
308  hops2 += distance2;
309  for (int k = 0 ; k < mach_coord_dim ; ++k) {
310  part_t distance =
311  ZOLTAN2_ABS(proc_coords[k][procId1] - proc_coords[k][procId2]);
312  if (machine_extent_wrap_around[k]) {
313  if (machine_extent[k] - distance < distance) {
314  distance = machine_extent[k] - distance;
315  }
316  }
317  hops += distance;
318  }
319  }
320  }
321  delete [] machine_extent_wrap_around;
322  delete [] machine_extent;
323 
324  if (tcomm->getRank() == 0)
325  std::cout << "HOPS:" << hops << " HOPS2:" << hops2 << std::endl;
326 
327  delete [] task_communication_xadj_tmp;
328  delete [] task_communication_adj_tmp;
329  }
330  /*
331  part_t *machineDimensions = NULL;
332  //machineDimensions = hopper;
333  Zoltan2::coordinateTaskMapperInterface<part_t, zscalar_t, zscalar_t>(
334  tcomm,
335  procDim,
336  numProcs,
337  procCoordinates,
338 
339  coordDim,
340  numGlobalTasks,
341  partCenters,
342 
343  task_communication_xadj_,
344  task_communication_adj_,
345  NULL,
346 
347  proc_to_task_xadj_,
348  proc_to_task_adj_,
349 
350  partArraysize,
351  partArray,
352  machineDimensions
353  );
354  */
355 
356  if (tcomm->getRank() == 0) {
357  std::cout << "PASS" << std::endl;
358  }
359 
360 
361  for (int i = 0; i < coordDim; i++) delete [] partCenters[i];
362  delete [] partCenters;
363 
364  }
365  catch(std::string &s) {
366  std::cerr << s << std::endl;
367  }
368 
369  catch(char * s) {
370  std::cerr << s << std::endl;
371  }
372 }
373 
Zoltan2::XpetraCrsGraphAdapter< tcrsGraph_t, tMVector_t > my_adapter_t
bool getHopCount(int rank1, int rank2, pcoord_t &hops) const
return the hop count between rank1 and rank2
int main(int narg, char *arg[])
double zscalar_t
Provides access for Zoltan2 to Xpetra::CrsGraph data.
int zlno_t
common code used by tests
Tpetra::Map::node_type mytest_znode_t
part_t getAssignedProcForTask(part_t taskId)
getAssignedProcForTask function, returns the assigned processor id for the given task ...
Defines the XpetraMultiVectorAdapter.
Defines XpetraCrsGraphAdapter class.
SparseMatrixAdapter_t::part_t part_t
int getMachineDim() const
returns the dimension (number of coords per node) in the machine
bool getAllMachineCoordinatesView(pcoord_t **&allCoords) const
getProcDim function set the coordinates of all ranks allCoords[i][j], i=0,...,getMachineDim(), j=0,...,getNumRanks(), is the i-th dimensional coordinate for rank j. return true if coordinates are available for all ranks
void getProcTask(part_t *&proc_to_task_xadj_, part_t *&proc_to_task_adj_)
#define ZOLTAN2_ABS(x)
Tpetra::CrsGraph< zlno_t, zgno_t, znode_t > tcrsGraph_t
Definition: GraphModel.cpp:83
bool getMachineExtentWrapArounds(bool *wrap_around) const
if the machine has a wrap-around tourus link in each dimension.
InputTraits< User >::part_t part_t
An adapter for Xpetra::MultiVector.
MachineRepresentation Class Base class for representing machine coordinates, networks, etc.
int zgno_t
PartitioningProblem sets up partitioning problems for the user.
Tpetra::MultiVector< double, int, int > tMVector_t
Defines the PartitioningProblem class.
A class that computes and returns quality metrics.
bool getMachineExtent(int *nxyz) const
sets the number of unique coordinates in each machine dimension