Zoltan2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
rcbPerformanceZ1.cpp
Go to the documentation of this file.
1 // @HEADER
2 //
3 // ***********************************************************************
4 //
5 // Zoltan2: A package of combinatorial algorithms for scientific computing
6 // Copyright 2012 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact Karen Devine (kddevin@sandia.gov)
39 // Erik Boman (egboman@sandia.gov)
40 // Siva Rajamanickam (srajama@sandia.gov)
41 //
42 // ***********************************************************************
43 //
44 // @HEADER
45 
53 #include <Zoltan2_TestHelpers.hpp>
54 
55 #include <zoltan.h>
56 #include <Teuchos_CommandLineProcessor.hpp>
57 
58 #include <vector>
59 #include <ostream>
60 #include <sstream>
61 #include <string>
62 
63 #include <Zoltan2_TestHelpers.hpp>
67 #include <GeometricGenerator.hpp>
69 
70 using std::vector;
71 using std::bad_alloc;
72 using Teuchos::RCP;
73 using Teuchos::rcp;
74 using Teuchos::Comm;
75 using Teuchos::ArrayView;
76 using Teuchos::ArrayRCP;
77 using Teuchos::CommandLineProcessor;
78 
79 typedef Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t> tMVector_t;
80 typedef Tpetra::Map<zlno_t, zgno_t, znode_t> tMap_t;
81 
82 static ArrayRCP<ArrayRCP<zscalar_t> > weights;
83 static RCP<tMVector_t> coordinates;
84 
85 
86 const char param_comment = '#';
87 
88 std::string trim_right_copy(
89  const std::string& s,
90  const std::string& delimiters = " \f\n\r\t\v" )
91 {
92  return s.substr( 0, s.find_last_not_of( delimiters ) + 1 );
93 }
94 
95 std::string trim_left_copy(
96  const std::string& s,
97  const std::string& delimiters = " \f\n\r\t\v" )
98 {
99  return s.substr( s.find_first_not_of( delimiters ) );
100 }
101 
102 std::string trim_copy(
103  const std::string& s,
104  const std::string& delimiters = " \f\n\r\t\v" )
105 {
106  return trim_left_copy( trim_right_copy( s, delimiters ), delimiters );
107 }
108 
109 // Zoltan1 query functions
110 bool getArgumentValue(std::string &argumentid, double &argumentValue, std::string argumentline){
111  std::stringstream stream(std::stringstream::in | std::stringstream::out);
112  stream << argumentline;
113  getline(stream, argumentid, '=');
114  if (stream.eof()){
115  return false;
116  }
117  stream >> argumentValue;
118  return true;
119 }
120 
121 std::string convert_to_string(char *args){
122  std::string tmp = "";
123  for(int i = 0; args[i] != 0; i++)
124  tmp += args[i];
125  return tmp;
126 }
127 int getNumObj(void *data, int *ierr)
128 {
129  *ierr = 0;
130 
131  return coordinates->getLocalLength();
132 }
133 
134 int getDim(void *data, int *ierr)
135 {
136  *ierr = 0;
137  return 3;
138 }
139 
140 void getObjList(void *data, int numGid, int numLid,
141  zgno_t * gids, zgno_t * lids,
142  int num_wgts, float *obj_wgts, int *ierr)
143 {
144  *ierr = 0;
145  size_t localLen = coordinates->getLocalLength();
146  const zgno_t *ids = coordinates->getMap()->getNodeElementList().getRawPtr();
147  zgno_t *idsNonConst = const_cast<zgno_t *>(ids);
148 
149  if (sizeof(zgno_t) == sizeof(zgno_t)){
150  memcpy(gids, idsNonConst, sizeof(zgno_t) * localLen);
151  }
152  else{
153  for (size_t i=0; i < localLen; i++)
154  gids[i] = static_cast<zgno_t>(idsNonConst[i]);
155  }
156 
157  if (num_wgts > 0){
158  float *wgts = obj_wgts;
159  for (size_t i=0; i < localLen; i++)
160  for (int w=0; w < num_wgts; w++)
161  *wgts++ = static_cast<float>(weights[w][i]);
162  }
163 }
164 
165 void getCoords(void *data, int numGid, int numLid,
166  int numObj, zgno_t * gids, zgno_t * lids,
167  int dim, double *coords, int *ierr)
168 {
169  // I know that Zoltan asks for coordinates in gid order.
170  *ierr = 0;
171  double *val = coords;
172  const zscalar_t *x = coordinates->getData(0).getRawPtr();
173  const zscalar_t *y = coordinates->getData(1).getRawPtr();
174  const zscalar_t *z = coordinates->getData(2).getRawPtr();
175  for (zlno_t i=0; i < numObj; i++){
176  *val++ = static_cast<double>(x[i]);
177  *val++ = static_cast<double>(y[i]);
178  *val++ = static_cast<double>(z[i]);
179  }
180 }
181 
182 
183 
189 };
190 
191 ArrayRCP<zscalar_t> makeWeights(
192  const RCP<const Teuchos::Comm<int> > & comm,
193  zlno_t len, weightTypes how, zscalar_t scale, int rank)
194 {
195  zscalar_t *wgts = new zscalar_t [len];
196  if (!wgts)
197  throw bad_alloc();
198 
199  ArrayRCP<zscalar_t> wgtArray(wgts, 0, len, true);
200 
201  if (how == upDown){
202  zscalar_t val = scale + rank%2;
203  for (zlno_t i=0; i < len; i++)
204  wgts[i] = val;
205  }
206  else if (how == roundRobin){
207  for (int i=0; i < 10; i++){
208  zscalar_t val = (i + 10)*scale;
209  for (int j=i; j < len; j += 10)
210  wgts[j] = val;
211  }
212  }
213  else if (how == increasing){
214  zscalar_t val = scale + rank;
215  for (zlno_t i=0; i < len; i++)
216  wgts[i] = val;
217  }
218 
219  return wgtArray;
220 }
221 
226 const RCP<tMVector_t> getMeshCoordinates(
227  const RCP<const Teuchos::Comm<int> > & comm,
228  zgno_t numGlobalCoords)
229 {
230  int rank = comm->getRank();
231  int nprocs = comm->getSize();
232 
233  double k = log(numGlobalCoords) / 3;
234  double xdimf = exp(k) + 0.5;
235  zgno_t xdim = static_cast<int>(floor(xdimf));
236  zgno_t ydim = xdim;
237  zgno_t zdim = numGlobalCoords / (xdim*ydim);
238  zgno_t num=xdim*ydim*zdim;
239  zgno_t diff = numGlobalCoords - num;
240  zgno_t newdiff = 0;
241 
242  while (diff > 0){
243  if (zdim > xdim && zdim > ydim){
244  zdim++;
245  newdiff = diff - (xdim*ydim);
246  if (newdiff < 0)
247  if (diff < -newdiff)
248  zdim--;
249  }
250  else if (ydim > xdim && ydim > zdim){
251  ydim++;
252  newdiff = diff - (xdim*zdim);
253  if (newdiff < 0)
254  if (diff < -newdiff)
255  ydim--;
256  }
257  else{
258  xdim++;
259  newdiff = diff - (ydim*zdim);
260  if (newdiff < 0)
261  if (diff < -newdiff)
262  xdim--;
263  }
264 
265  diff = newdiff;
266  }
267 
268  num=xdim*ydim*zdim;
269  diff = numGlobalCoords - num;
270  if (diff < 0)
271  diff /= -numGlobalCoords;
272  else
273  diff /= numGlobalCoords;
274 
275  if (rank == 0){
276  if (diff > .01)
277  std::cout << "Warning: Difference " << diff*100 << " percent" << std::endl;
278  std::cout << "Mesh size: " << xdim << "x" << ydim << "x" <<
279  zdim << ", " << num << " vertices." << std::endl;
280  }
281 
282  // Divide coordinates.
283 
284  zgno_t numLocalCoords = num / nprocs;
285  zgno_t leftOver = num % nprocs;
286  zgno_t gid0 = 0;
287 
288  if (rank <= leftOver)
289  gid0 = zgno_t(rank) * (numLocalCoords+1);
290  else
291  gid0 = (leftOver * (numLocalCoords+1)) +
292  ((zgno_t(rank) - leftOver) * numLocalCoords);
293 
294  if (rank < leftOver)
295  numLocalCoords++;
296 
297  zgno_t gid1 = gid0 + numLocalCoords;
298 
299  zgno_t *ids = new zgno_t [numLocalCoords];
300  if (!ids)
301  throw bad_alloc();
302  ArrayRCP<zgno_t> idArray(ids, 0, numLocalCoords, true);
303 
304  for (zgno_t i=gid0; i < gid1; i++)
305  *ids++ = i;
306 
307  RCP<const tMap_t> idMap = rcp(
308  new tMap_t(num, idArray.view(0, numLocalCoords), 0, comm));
309 
310  // Create a Tpetra::MultiVector of coordinates.
311 
312  zscalar_t *x = new zscalar_t [numLocalCoords*3];
313  if (!x)
314  throw bad_alloc();
315  ArrayRCP<zscalar_t> coordArray(x, 0, numLocalCoords*3, true);
316 
317  zscalar_t *y = x + numLocalCoords;
318  zscalar_t *z = y + numLocalCoords;
319 
320  zgno_t xStart = 0;
321  zgno_t yStart = 0;
322  zgno_t xyPlane = xdim*ydim;
323  zgno_t zStart = gid0 / xyPlane;
324  zgno_t rem = gid0 % xyPlane;
325  if (rem > 0){
326  yStart = rem / xdim;
327  xStart = rem % xdim;
328  }
329 
330  zlno_t next = 0;
331  for (zscalar_t zval=zStart; next < numLocalCoords && zval < zdim; zval++){
332  for (zscalar_t yval=yStart; next < numLocalCoords && yval < ydim; yval++){
333  for (zscalar_t xval=xStart; next < numLocalCoords && xval < xdim; xval++){
334  x[next] = xval;
335  y[next] = yval;
336  z[next] = zval;
337  next++;
338  }
339  xStart = 0;
340  }
341  yStart = 0;
342  }
343 
344  ArrayView<const zscalar_t> xArray(x, numLocalCoords);
345  ArrayView<const zscalar_t> yArray(y, numLocalCoords);
346  ArrayView<const zscalar_t> zArray(z, numLocalCoords);
347  ArrayRCP<ArrayView<const zscalar_t> > coordinates =
348  arcp(new ArrayView<const zscalar_t> [3], 0, 3);
349  coordinates[0] = xArray;
350  coordinates[1] = yArray;
351  coordinates[2] = zArray;
352 
353  ArrayRCP<const ArrayView<const zscalar_t> > constCoords =
354  coordinates.getConst();
355 
356  RCP<tMVector_t> meshCoords = rcp(new tMVector_t(
357  idMap, constCoords.view(0,3), 3));
358 
359  return meshCoords;
360 }
361 
362 
363 void getArgVals(int narg, char **arg, int &numParts,
364  std::string &paramFile){
365 
366  for(int i = 0; i < narg; ++i){
367  std::string tmp = convert_to_string(arg[i]);
368  std::string identifier = "";
369  long long int value = -1; double fval = -1;
370  if(!getArgumentValue(identifier, fval, tmp)) continue;
371  value = (long long int) (fval);
372  if(identifier == "C"){
373  if(value > 0){
374  numParts=value;
375  } else {
376  throw "Invalid argument at " + tmp;
377  }
378  }
379  else if(identifier == "PF"){
380  std::stringstream stream(std::stringstream::in | std::stringstream::out);
381  stream << tmp;
382  getline(stream, paramFile, '=');
383  stream >> paramFile;
384  }
385  else {
386  throw "Invalid argument at " + tmp;
387  }
388 
389  }
390 
391 }
392 
393 void readGeoGenParams(std::string paramFileName, Teuchos::ParameterList &geoparams, const RCP<const Teuchos::Comm<int> > & comm){
394  std::string input = "";
395  char inp[25000];
396  for(int i = 0; i < 25000; ++i){
397  inp[i] = 0;
398  }
399  if(comm->getRank() == 0){
400  std::fstream inParam(paramFileName.c_str());
401 
402  std::string tmp = "";
403  getline (inParam,tmp);
404  while (!inParam.eof()){
405  if(tmp != ""){
406  tmp = trim_copy(tmp);
407  if(tmp != ""){
408  input += tmp + "\n";
409  }
410  }
411  getline (inParam,tmp);
412  }
413  inParam.close();
414  for (size_t i = 0; i < input.size(); ++i){
415  inp[i] = input[i];
416  }
417  }
418 
419 
420  int size = input.size();
421 
422  //MPI_Bcast(&size,1,MPI_INT, 0, MPI_COMM_WORLD);
423  //MPI_Bcast(inp,size,MPI_CHAR, 0, MPI_COMM_WORLD);
424  //Teuchos::broadcast<int, char>(comm, 0,inp);
425 
426  comm->broadcast(0, sizeof(int), (char*) &size);
427  comm->broadcast(0, size, inp);
428  //Teuchos::broadcast<int,std::string>(comm,0, &input);
429  std::istringstream inParam(inp);
430  std::string str;
431  getline (inParam,str);
432  while (!inParam.eof()){
433  if(str[0] != param_comment){
434  size_t pos = str.find('=');
435  if(pos == std::string::npos){
436  throw "Invalid Line:" + str + " in parameter file";
437  }
438  std::string paramname = trim_copy(str.substr(0,pos));
439  std::string paramvalue = trim_copy(str.substr(pos + 1));
440  geoparams.set(paramname, paramvalue);
441  }
442  getline (inParam,str);
443  }
444 }
445 
446 int main(int narg, char *arg[])
447 {
448  // MEMORY_CHECK(true, "Before initializing MPI");
449  Tpetra::ScopeGuard tscope(&narg, &arg);
450  Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
451 
452  int rank = comm->getRank();
453  int nprocs = comm->getSize();
454 
455  MEMORY_CHECK(rank==0 || rank==nprocs-1, "After initializing MPI");
456 
457  if (rank==0)
458  std::cout << "Number of processes: " << nprocs << std::endl;
459 
460  // Default values
461  //double numGlobalCoords = 1000;
462  //int nWeights = 0;
463  int debugLevel=2; // for timing
464  std::string memoryOn("memoryOn");
465  std::string memoryOff("memoryOff");
466  bool doMemory=false;
467 
468  int numGlobalParts = 512 * 512;
469  std::string paramFile = "p2.txt";
470  //getArgVals(narg, arg, numGlobalParts, paramFile);
471 
472  int dummyTimer=0;
473 
474 
475 
476  Teuchos::ParameterList geoparams("geo params");
477 
478  readGeoGenParams(paramFile, geoparams, comm);
479 #ifdef HAVE_ZOLTAN2_OMP
480  double begin = omp_get_wtime();
481 #endif
482  GeometricGenerator<zscalar_t, zlno_t, zgno_t, znode_t> *gg = new GeometricGenerator<zscalar_t, zlno_t, zgno_t, znode_t>(geoparams,comm);
483 #ifdef HAVE_ZOLTAN2_OMP
484  double end = omp_get_wtime();
485  std::cout << "GeometricGen Time:" << end - begin << std::endl;
486 #endif
487  int coord_dim = gg->getCoordinateDimension();
488  int nWeights = gg->getNumWeights();
489  zlno_t numLocalPoints = gg->getNumLocalCoords();
490  zgno_t numGlobalPoints = gg->getNumGlobalCoords();
491  zscalar_t **coords = new zscalar_t * [coord_dim];
492  for(int i = 0; i < coord_dim; ++i){
493  coords[i] = new zscalar_t[numLocalPoints];
494  }
495  gg->getLocalCoordinatesCopy(coords);
496  zscalar_t **weight = NULL;
497  if(nWeights){
498  weight= new zscalar_t * [nWeights];
499  for(int i = 0; i < nWeights; ++i){
500  weight[i] = new zscalar_t[numLocalPoints];
501  }
502  gg->getLocalWeightsCopy(weight);
503  }
504 
505  zgno_t globalSize = static_cast<zgno_t>(numGlobalPoints);
506  delete gg;
507 
508  std::cout << "coord_dim:" << coord_dim << " nWeights:" << nWeights << " numLocalPoints:" << numLocalPoints << " numGlobalPoints:" << numGlobalPoints << std::endl;
509 
510  CommandLineProcessor commandLine(false, true);
511  commandLine.setOption("size", &numGlobalPoints,
512  "Approximate number of global coordinates.");
513  commandLine.setOption("numParts", &numGlobalParts,
514  "Number of parts (default is one per proc).");
515  commandLine.setOption("numWeights", &nWeights,
516  "Number of weights per coordinate, zero implies uniform weights.");
517  commandLine.setOption("debug", &debugLevel, "Zoltan1 debug level");
518  commandLine.setOption("timers", &dummyTimer, "ignored");
519  commandLine.setOption(memoryOn.c_str(), memoryOff.c_str(), &doMemory,
520  "do memory profiling");
521 
522  std::string balanceCount("balance_object_count");
523  std::string balanceWeight("balance_object_weight");
524  std::string mcnorm1("multicriteria_minimize_total_weight");
525  std::string mcnorm2("multicriteria_balance_total_maximum");
526  std::string mcnorm3("multicriteria_minimize_maximum_weight");
527 
528  std::string objective(balanceWeight); // default
529 
530  std::string doc(balanceCount);
531  doc.append(": ignore weights\n");
532 
533  doc.append(balanceWeight);
534  doc.append(": balance on first weight\n");
535 
536  doc.append(mcnorm1);
537  doc.append(": given multiple weights, balance their total.\n");
538 
539  doc.append(mcnorm3);
540  doc.append(": given multiple weights, balance the maximum for each coordinate.\n");
541 
542  doc.append(mcnorm2);
543  doc.append(": given multiple weights, balance the L2 norm of the weights.\n");
544 
545  commandLine.setOption("objective", &objective, doc.c_str());
546 
547  CommandLineProcessor::EParseCommandLineReturn rc =
548  commandLine.parse(narg, arg);
549 
550  if (rc != Teuchos::CommandLineProcessor::PARSE_SUCCESSFUL){
551  if (rc == Teuchos::CommandLineProcessor::PARSE_HELP_PRINTED){
552  if (rank==0)
553  std::cout << "PASS" << std::endl;
554  return 1;
555  }
556  else{
557  if (rank==0)
558  std::cout << "FAIL" << std::endl;
559  return 0;
560  }
561  }
562 
563  //MEMORY_CHECK(doMemory && rank==0, "After processing parameters");
564 
565 
566 
567 
568  RCP<Tpetra::Map<zlno_t, zgno_t, znode_t> > mp = rcp(
569  new Tpetra::Map<zlno_t, zgno_t, znode_t> (numGlobalPoints, numLocalPoints, 0, comm));
570 
571  Teuchos::Array<Teuchos::ArrayView<const zscalar_t> > coordView(coord_dim);
572  for (int i=0; i < coord_dim; i++){
573  if(numLocalPoints > 0){
574  Teuchos::ArrayView<const zscalar_t> a(coords[i], numLocalPoints);
575  coordView[i] = a;
576  } else{
577  Teuchos::ArrayView<const zscalar_t> a;
578  coordView[i] = a;
579  }
580  }
581 
582  coordinates = RCP< Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t> >(
583  new Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t>( mp, coordView.view(0, coord_dim), coord_dim));
584 
585  //typedef Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t> tMVector_t;
586 
587  RCP<const tMVector_t> coordsConst = Teuchos::rcp_const_cast<const tMVector_t>(coordinates);
588 
589 
590 
591  //coordinates = getMeshCoordinates(comm, globalSize);
592  size_t numLocalCoords = coordinates->getLocalLength();
593 
594 #if 0
595  comm->barrier();
596  for (int p=0; p < nprocs; p++){
597  if (p==rank){
598  std::cout << "Rank " << rank << ", " << numLocalCoords << "coords" << std::endl;
599  const zscalar_t *x = coordinates->getData(0).getRawPtr();
600  const zscalar_t *y = coordinates->getData(1).getRawPtr();
601  const zscalar_t *z = coordinates->getData(2).getRawPtr();
602  for (zlno_t i=0; i < numLocalCoords; i++)
603  std::cout << " " << x[i] << " " << y[i] << " " << z[i] << std::endl;
604  }
605  std::cout.flush();
606  comm->barrier();
607  }
608 #endif
609 
610  if (nWeights > 0){
611 
612  weights = arcp(new ArrayRCP<zscalar_t> [nWeights],
613  0, nWeights, true);
614  for(int i = 0; i < nWeights; ++i){
615  //weights[i] = ArrayRCP<zscalar_t>(weight[i]);
616  }
617  }
618 
619  MEMORY_CHECK(doMemory && rank==0, "After creating input");
620 
621  // Now call Zoltan to partition the problem.
622 
623  float ver;
624  int aok = Zoltan_Initialize(narg, arg, &ver);
625 
626  if (aok != 0){
627  printf("sorry...\n");
628  exit(0);
629  }
630 
631  struct Zoltan_Struct *zz;
632  zz = Zoltan_Create(MPI_COMM_WORLD);
633 
634  Zoltan_Set_Param(zz, "LB_METHOD", "RCB");
635  Zoltan_Set_Param(zz, "LB_APPROACH", "PARTITION");
636  Zoltan_Set_Param(zz, "CHECK_GEOM", "0");
637  Zoltan_Set_Param(zz, "NUM_GID_ENTRIES", "1"); // compiled with ULONG option
638  Zoltan_Set_Param(zz, "NUM_LID_ENTRIES", "0");
639  Zoltan_Set_Param(zz, "RETURN_LISTS", "ALL");
640 
641  //Zoltan_Set_Param(zz, "FINAL_OUTPUT", "1");
642  Zoltan_Set_Param(zz, "IMBALANCE_TOL", "1.0");
643  std::ostringstream oss;
644  oss << numGlobalParts;
645  Zoltan_Set_Param(zz, "NUM_GLOBAL_PARTS", oss.str().c_str());
646  oss.str("");
647  oss << debugLevel;
648  Zoltan_Set_Param(zz, "DEBUG_LEVEL", oss.str().c_str());
649 
650  if (objective != balanceCount){
651  oss.str("");
652  oss << nWeights;
653  Zoltan_Set_Param(zz, "OBJ_WEIGHT_DIM", oss.str().c_str());
654 
655  if (objective == mcnorm1)
656  Zoltan_Set_Param(zz, "RCB_MULTICRITERIA_NORM", "1");
657  else if (objective == mcnorm2)
658  Zoltan_Set_Param(zz, "RCB_MULTICRITERIA_NORM", "2");
659  else if (objective == mcnorm3)
660  Zoltan_Set_Param(zz, "RCB_MULTICRITERIA_NORM", "3");
661  }
662  else{
663  Zoltan_Set_Param(zz, "OBJ_WEIGHT_DIM", "0");
664  }
665 
666  Zoltan_Set_Num_Obj_Fn(zz, getNumObj, NULL);
667  Zoltan_Set_Obj_List_Fn(zz, getObjList,NULL);
668  Zoltan_Set_Num_Geom_Fn(zz, getDim, NULL);
669  Zoltan_Set_Geom_Multi_Fn(zz, getCoords, NULL);
670 
671  int changes, numGidEntries, numLidEntries, numImport, numExport;
672  zgno_t * importGlobalGids, * importLocalGids;
673  zgno_t * exportGlobalGids, * exportLocalGids;
674  int *importProcs, *importToPart, *exportProcs, *exportToPart;
675 
676  MEMORY_CHECK(doMemory && rank==0, "Before Zoltan_LB_Partition");
677 
678 
679  aok = Zoltan_LB_Partition(zz, /* input (all remaining fields are output) */
680  &changes, /* 1 if partitioning was changed, 0 otherwise */
681  &numGidEntries, /* Number of integers used for a global ID */
682  &numLidEntries, /* Number of integers used for a local ID */
683  &numImport, /* Number of vertices to be sent to me */
684  &importGlobalGids, /* Global IDs of vertices to be sent to me */
685  &importLocalGids, /* Local IDs of vertices to be sent to me */
686  &importProcs, /* Process rank for source of each incoming vertex */
687  &importToPart, /* New partition for each incoming vertex */
688  &numExport, /* Number of vertices I must send to other processes*/
689  &exportGlobalGids, /* Global IDs of the vertices I must send */
690  &exportLocalGids, /* Local IDs of the vertices I must send */
691  &exportProcs, /* Process to which I send each of the vertices */
692  &exportToPart); /* Partition to which each vertex will belong */
693 
694  MEMORY_CHECK(doMemory && rank==0, "After Zoltan_LB_Partition");
695  Zoltan_LB_Free_Part(importGlobalGids, importLocalGids, importProcs, importToPart);
696  Zoltan_LB_Free_Part(exportGlobalGids, exportLocalGids, exportProcs, exportToPart);
697  Zoltan_Destroy(&zz);
698  MEMORY_CHECK(doMemory && rank==0, "After Zoltan_Destroy");
699 
700  if (rank==0){
701  if (aok != 0)
702  std::cout << "FAIL" << std::endl;
703  else
704  std::cout << "PASS" << std::endl;
705  }
706 
707  return 0;
708 }
string trim_right_copy(const string &s, const string &delimiters=" \f\n\r\t\v")
Tpetra::Map< zlno_t, zgno_t, znode_t > tMap_t
void readGeoGenParams(string paramFileName, Teuchos::ParameterList &geoparams, const RCP< const Teuchos::Comm< int > > &comm)
int main(int narg, char *arg[])
double zscalar_t
static ArrayRCP< ArrayRCP< zscalar_t > > weights
int zlno_t
Defines the PartitioningSolution class.
common code used by tests
int getNumObj(void *data, int *ierr)
static RCP< tMVector_t > coordinates
Defines the XpetraMultiVectorAdapter.
string trim_copy(const string &s, const string &delimiters=" \f\n\r\t\v")
const char param_comment
void getCoords(void *data, int numGid, int numLid, int numObj, zgno_t *gids, zgno_t *lids, int dim, double *coords, int *ierr)
Defines the EvaluatePartition class.
int getDim(void *data, int *ierr)
int zgno_t
string convert_to_string(char *args)
ArrayRCP< zscalar_t > makeWeights(const RCP< const Teuchos::Comm< int > > &comm, zlno_t len, weightTypes how, zscalar_t scale, int rank)
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)
Defines the PartitioningProblem class.
void getObjList(void *data, int numGid, int numLid, zgno_t *gids, zgno_t *lids, int num_wgts, float *obj_wgts, int *ierr)
const RCP< tMVector_t > getMeshCoordinates(const RCP< const Teuchos::Comm< int > > &comm, zgno_t numGlobalCoords)
Create a mesh of approximately the desired size.
#define MEMORY_CHECK(iPrint, msg)