Zoltan2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
TaskMappingTest3.cpp
Go to the documentation of this file.
3 #include "Tpetra_MultiVector_decl.hpp"
4 
5 #include "GeometricGenerator.hpp"
6 #include <string>
7 #include "Teuchos_XMLParameterListHelpers.hpp"
8 //#include "Teuchos_MPIComm.hpp"
9 #define partDIM 3
10 #define procDIM 3
11 #define nProcs 200;
12 #define nParts 200;
13 
14 
16  const string& s,
17  const string& delimiters = " \f\n\r\t\v" )
18 {
19  return s.substr( 0, s.find_last_not_of( delimiters ) + 1 );
20 }
21 
23  const string& s,
24  const string& delimiters = " \f\n\r\t\v" )
25 {
26  return s.substr( s.find_first_not_of( delimiters ) );
27 }
28 
29 string trim_copy(
30  const string& s,
31  const string& delimiters = " \f\n\r\t\v" )
32 {
33  return trim_left_copy( trim_right_copy( s, delimiters ), delimiters );
34 }
35 
36 const char param_comment = '#';
37 void readGeoGenParams(string paramFileName, Teuchos::ParameterList &geoparams, const RCP<const Teuchos::Comm<int> > & comm){
38  std::string input = "";
39  char inp[25000];
40  for(int i = 0; i < 25000; ++i){
41  inp[i] = 0;
42  }
43 
44  bool fail = false;
45  if(comm->getRank() == 0){
46 
47  std::fstream inParam(paramFileName.c_str());
48  if (inParam.fail())
49  {
50  fail = true;
51  }
52  if(!fail)
53  {
54  std::string tmp = "";
55  getline (inParam,tmp);
56  while (!inParam.eof()){
57  if(tmp != ""){
58  tmp = trim_copy(tmp);
59  if(tmp != ""){
60  input += tmp + "\n";
61  }
62  }
63  getline (inParam,tmp);
64  }
65  inParam.close();
66  for (size_t i = 0; i < input.size(); ++i){
67  inp[i] = input[i];
68  }
69  }
70  }
71 
72 
73 
74  int size = input.size();
75  if(fail){
76  size = -1;
77  }
78  comm->broadcast(0, sizeof(int), (char*) &size);
79  if(size == -1){
80  throw "File " + paramFileName + " cannot be opened.";
81  }
82  comm->broadcast(0, size, inp);
83  std::istringstream inParam(inp);
84  std::string str;
85  getline (inParam,str);
86  while (!inParam.eof()){
87  if(str[0] != param_comment){
88  size_t pos = str.find('=');
89  if(pos == std::string::npos){
90  throw "Invalid Line:" + str + " in parameter file";
91  }
92  string paramname = trim_copy(str.substr(0,pos));
93  string paramvalue = trim_copy(str.substr(pos + 1));
94  geoparams.set(paramname, paramvalue);
95  }
96  getline (inParam,str);
97  }
98 }
99 string convert_to_string(char *args){
100  string tmp = "";
101  for(int i = 0; args[i] != 0; i++)
102  tmp += args[i];
103  return tmp;
104 }
105 bool getArgumentValue(string &argumentid, double &argumentValue, string argumentline){
106  std::stringstream stream(std::stringstream::in | std::stringstream::out);
107  stream << argumentline;
108  getline(stream, argumentid, '=');
109  if (stream.eof()){
110  return false;
111  }
112  stream >> argumentValue;
113  return true;
114 }
115 
116 template <typename part_t>
118  int narg,
119  char **argv,
120  std::string &procF,
121  part_t &nx,
122  part_t &ny,
123  part_t &nz, bool &divide_prime, int &ranks_per_node, string &taskGraphFile, string &taskCoordFile){
124 
125  bool isprocset = false;
126  int ispartset = 0;
127 
128  for(int i = 0; i < narg; ++i){
129  string tmp = convert_to_string(argv[i]);
130  string tmp2 = "";
131  string identifier = "";
132  double fval = -1;
133  if(!getArgumentValue(identifier, fval, tmp)) continue;
134 
135  if(identifier == "PROC"){
136  std::stringstream stream(std::stringstream::in | std::stringstream::out);
137  stream << tmp;
138  getline(stream, procF, '=');
139 
140  stream >> procF;
141  isprocset = true;
142  }
143  else if(identifier == "NX"){
144  std::stringstream stream(std::stringstream::in | std::stringstream::out);
145  stream << tmp;
146  getline(stream, tmp2, '=');
147 
148  stream >> nx;
149  ispartset++;
150  }
151  else if(identifier == "NY"){
152  std::stringstream stream(std::stringstream::in | std::stringstream::out);
153  stream << tmp;
154  getline(stream, tmp2, '=');
155 
156  stream >> ny;
157  ispartset++;
158  }
159  else if(identifier == "NZ"){
160  std::stringstream stream(std::stringstream::in | std::stringstream::out);
161  stream << tmp;
162  getline(stream, tmp2, '=');
163 
164  stream >> nz;
165  ispartset++;
166  }
167  else if(identifier == "DP"){
168  std::stringstream stream(std::stringstream::in | std::stringstream::out);
169  stream << tmp;
170  getline(stream, tmp2, '=');
171  int val;
172  stream >> val;
173  if (val) divide_prime = true;
174  ispartset++;
175  }
176  else if(identifier == "RPN"){
177  std::stringstream stream(std::stringstream::in | std::stringstream::out);
178  stream << tmp;
179  getline(stream, tmp2, '=');
180  stream >> ranks_per_node;
181  ispartset++;
182  } else if(identifier == "TG"){
183  std::stringstream stream(std::stringstream::in | std::stringstream::out);
184  stream << tmp;
185  getline(stream, taskGraphFile, '=');
186 
187  stream >> taskGraphFile;
188  }
189  else if(identifier == "TC"){
190  std::stringstream stream(std::stringstream::in | std::stringstream::out);
191  stream << tmp;
192  getline(stream, taskCoordFile, '=');
193  stream >> taskCoordFile;
194  }
195 
196  else {
197  throw "Invalid argument at " + tmp;
198  }
199 
200  }
201  if(!(ispartset >= 3&& isprocset)){
202  throw "(PROC && PART) are mandatory arguments.";
203  }
204 
205 }
206 int main(int narg, char *arg[]){
207 
208  Tpetra::ScopeGuard tscope(&narg, &arg);
209 
210  typedef Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t> tMVector_t;
211  typedef Zoltan2::XpetraMultiVectorAdapter<tMVector_t> inputAdapter_t;
213 
214  //if (narg != 3){
215  // std::cout << "Usage: " << arg[0] << " PART=partGeoParams.txt PROC=procGeoParams.txt" << std::endl;
216  // exit(1);
217  //}
218  part_t numParts = 0;
219  zscalar_t **partCenters = NULL;
220  int coordDim = 0;
221 
222  part_t numProcs = 0;
223  zscalar_t **procCoordinates = NULL;
224  int procDim = 0;
225 
226 
227  string taskGraphFile = "";
228  string taskCoordFile = "";
229  bool divide_prime = false;
230  part_t jobX = 1, jobY = 1, jobZ = 1;
231  string procfile = "";
232  int rank_per_node = 16;
233 
234  const RCP<Comm<int> > commN;
235  RCP<Comm<int> >comm = Teuchos::rcp_const_cast<Comm<int> >
236  (Teuchos::DefaultComm<int>::getDefaultSerialComm(commN));
237 
238  part_t *task_communication_xadj_ = NULL;
239  part_t *task_communication_adj_ = NULL;
240  zscalar_t *task_communication_adjw_ = NULL;
241  try {
242 
243  getArgVals<part_t>(
244  narg,
245  arg,
246  procfile ,
247  jobX, jobY, jobZ, divide_prime, rank_per_node, taskGraphFile, taskCoordFile);
248 
249  coordDim = 3;
250  procDim = 3;
251  numParts = jobZ*jobY*jobX;
252 
253  //numProcs = numParts;
254  //std::cout << "part:" << numParts << " proc:" << procfile << std::endl;
255  if (taskGraphFile == "" || taskCoordFile == "")
256  {
257 
258  partCenters = new zscalar_t * [coordDim];
259  for(int i = 0; i < coordDim; ++i){
260  partCenters[i] = new zscalar_t[numParts];
261  }
262 
263 
264  task_communication_xadj_ = new part_t [numParts+1];
265  task_communication_adj_ = new part_t [numParts * 6];
266 
267  int prevNCount = 0;
268  task_communication_xadj_[0] = 0;
269  for (part_t i = 0; i < numParts; ++i) {
270  int x = i % jobX;
271  int y = (i / (jobX)) % jobY;
272  int z = (i / (jobX)) / jobY;
273  partCenters[0][i] = x;
274  partCenters[1][i] = y;
275  partCenters[2][i] = z;
276 
277  if (x > 0){
278  task_communication_adj_[prevNCount++] = i - 1;
279  }
280  if (x < jobX - 1){
281  task_communication_adj_[prevNCount++] = i + 1;
282  }
283  if (y > 0){
284  task_communication_adj_[prevNCount++] = i - jobX;
285  }
286  if (y < jobY - 1){
287  task_communication_adj_[prevNCount++] = i + jobX;
288  }
289  if (z > 0){
290  task_communication_adj_[prevNCount++] = i - jobX * jobY;
291  }
292  if (z < jobZ - 1){
293  task_communication_adj_[prevNCount++] = i + jobX * jobY;
294  }
295  task_communication_xadj_[i+1] = prevNCount;
296  }
297  }
298  else {
299  int ne = 0;
300  FILE *f2 = fopen(taskGraphFile.c_str(), "rb");
301 
302  fread(&numParts,sizeof(int),1,f2); // write 10 bytes to our buffer
303  fread(&ne,sizeof(int),1,f2); // write 10 bytes to our buffer
304 
305 
306  std::cout << "numParts:" << numParts << " ne:" << ne << std::endl;
307 
308  task_communication_xadj_ = new part_t [numParts+1];
309  task_communication_adj_ = new part_t [ne];
310  task_communication_adjw_ = new zscalar_t [ne];
311 
312  fread((void *)task_communication_xadj_,sizeof(int),numParts + 1,f2); // write 10 bytes to our buffer
313  fread((void *)task_communication_adj_,sizeof(int),ne ,f2); // write 10 bytes to our buffer
314  fread((void *)task_communication_adjw_,sizeof(double),ne,f2); // write 10 bytes to our buffer
315  fclose(f2);
316 
317  f2 = fopen(taskCoordFile.c_str(), "rb");
318  fread((void *)&numParts,sizeof(int),1,f2); // write 10 bytes to our buffer
319  fread((void *)&coordDim,sizeof(int),1,f2); // write 10 bytes to our buffer
320 
321  std::cout << "numParts:" << numParts << " coordDim:" << coordDim << std::endl;
322 
323  partCenters = new zscalar_t * [coordDim];
324  for(int i = 0; i < coordDim; ++i){
325  partCenters[i] = new zscalar_t[numParts];
326  fread((void *) partCenters[i],sizeof(double),numParts, f2); // write 10 bytes to our buffer
327  }
328  fclose(f2);
329  }
330 
331 
332 
333 
334  {
335  std::vector < std::vector <zscalar_t> > proc_coords(procDim);
336  std::fstream m(procfile.c_str());
337  procCoordinates = new zscalar_t * [procDim];
338 
339  part_t i = 0;
340  zscalar_t a,b, c;
341  m >> a >> b >> c;
342  while(!m.eof()){
343  proc_coords[0].push_back(a);
344  proc_coords[1].push_back(b);
345  proc_coords[2].push_back(c);
346  ++i;
347  m >> a >> b >> c;
348  }
349 
350  m.close();
351  numProcs = i;
352  for(int ii = 0; ii < procDim; ++ii){
353  procCoordinates[ii] = new zscalar_t[numProcs];
354  for (int j = 0; j < numProcs; ++j){
355  procCoordinates[ii][j] = proc_coords[ii][j];
356  }
357  }
358  }
359  std::cout << "numProcs:" << numProcs << std::endl;
360 
361  /*
362  Zoltan2::CoordinateCommunicationModel<zscalar_t,zscalar_t,int> *cm =
363  new Zoltan2::CoordinateCommunicationModel<zscalar_t,zscalar_t,int>(
364  procDim, procCoordinates,
365  coordDim, partCenters,
366  numProcs, numParts);
367 
368  Zoltan2::Environment *env = new Zoltan2::Environment();
369  Zoltan2::CoordinateTaskMapper <inputAdapter_t, int> *ctm=
370  new Zoltan2::CoordinateTaskMapper<inputAdapter_t,int>(env, cm);
371 
372  */
373  Teuchos::RCP<const Teuchos::Comm<int> > tcomm =Tpetra::getDefaultComm();
374  part_t *proc_to_task_xadj_ = new part_t[numProcs+1];
375  part_t *proc_to_task_adj_ = new part_t[numParts];
376 /*
377  std::cout << "procDim:" << procDim <<
378  " numProcs:" << numProcs <<
379  " coordDim:" << coordDim <<
380  " numParts" << numParts << std::endl;
381 
382  for(part_t j = 0; j < numProcs; ++j){
383  std::cout << "proc - coord:" << j << " " << procCoordinates[0][j]<< " " << procCoordinates[1][j]<< " " << procCoordinates[2][j] << std::endl;
384  }
385 
386  for(part_t j = 0; j < numParts; ++j){
387  std::cout << "part - coord:" << j << " " << partCenters[0][j]<< " " << partCenters[1][j]<< " " << partCenters[2][j] << std::endl;
388  }
389 */
390  /*
391  int partArray[3];
392  partArray[0] = 8;
393  partArray[1] = 4;
394  partArray[2] = 16;
395  */
396  part_t *partArray = NULL;
397  int partArraysize = -1;
398  //part_t hopper[3];
399  //hopper[0] = 17;
400  //hopper[1] = 8;
401  //hopper[2] = 24;
402  part_t *machineDimensions = NULL;
403  //machineDimensions = hopper;
404 
405  Zoltan2::coordinateTaskMapperInterface<part_t, zscalar_t, zscalar_t>(
406  tcomm,
407  procDim,
408  numProcs,
409  procCoordinates,
410 
411  coordDim,
412  numParts,
413  partCenters,
414 
415  task_communication_xadj_,
416  task_communication_adj_,
417  task_communication_adjw_,
418 
419  proc_to_task_xadj_, /*output*/
420  proc_to_task_adj_, /*output*/
421 
422  partArraysize,
423  partArray,
424  machineDimensions, rank_per_node, divide_prime
425  );
426 
427  if (tcomm->getRank() == 0){
428  std::cout << "PASS" << std::endl;
429  }
430  /*
431  delete ctm;
432  delete cm;
433  delete env;
434  */
435  delete [] proc_to_task_xadj_;
436  delete [] proc_to_task_adj_;
437  delete [] task_communication_xadj_;
438  delete [] task_communication_adj_;
439  delete [] task_communication_adjw_;
440 
441  for (int i = 0; i < coordDim; i++) delete [] partCenters[i];
442  delete [] partCenters;
443  for (int i = 0; i < procDim; i++) delete [] procCoordinates[i];
444  delete [] procCoordinates;
445  }
446  catch(std::string &s){
447  std::cerr << s << std::endl;
448  }
449 
450  catch(char * s){
451  std::cerr << s << std::endl;
452  }
453 }
454 
string trim_right_copy(const string &s, const string &delimiters=" \f\n\r\t\v")
void readGeoGenParams(string paramFileName, Teuchos::ParameterList &geoparams, const RCP< const Teuchos::Comm< int > > &comm)
int main(int narg, char *arg[])
double zscalar_t
common code used by tests
string trim_copy(const string &s, const string &delimiters=" \f\n\r\t\v")
const char param_comment
SparseMatrixAdapter_t::part_t part_t
An adapter for Xpetra::MultiVector.
static const std::string fail
string convert_to_string(char *args)
string trim_left_copy(const string &s, const string &delimiters=" \f\n\r\t\v")
void getArgVals(int narg, char **arg, int &numParts, float &imbalance, string &pqParts, int &opt, std::string &fname, std::string &pfname, int &k, int &migration_check_option, int &migration_all_to_all_type, zscalar_t &migration_imbalance_cut_off, int &migration_processor_assignment_type, int &migration_doMigration_type, bool &test_boxes, bool &rectilinear, int &mj_premigration_option, int &mj_coordinate_cutoff)
Tpetra::MultiVector< double, int, int > tMVector_t
bool getArgumentValue(string &argumentid, double &argumentValue, string argumentline)