Zoltan2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
rcbPerformance.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 
52 #include <Zoltan2_TestHelpers.hpp>
56 
57 #include <Teuchos_CommandLineProcessor.hpp>
58 
59 #include <vector>
60 
61 using std::vector;
62 using std::bad_alloc;
63 using Teuchos::RCP;
64 using Teuchos::rcp;
65 using Teuchos::Comm;
66 using Teuchos::ArrayView;
67 using Teuchos::ArrayRCP;
68 using Teuchos::CommandLineProcessor;
69 
70 typedef Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t> tMVector_t;
71 typedef Tpetra::Map<zlno_t, zgno_t, znode_t> tMap_t;
72 
78 };
79 
80 ArrayRCP<zscalar_t> makeWeights(
81  const RCP<const Teuchos::Comm<int> > & comm,
82  zlno_t len, weightTypes how, zscalar_t scale, int rank)
83 {
84  zscalar_t *wgts = new zscalar_t [len];
85  if (!wgts)
86  throw bad_alloc();
87 
88  ArrayRCP<zscalar_t> weights(wgts, 0, len, true);
89 
90  if (how == upDown){
91  zscalar_t val = scale + rank%2;
92  for (zlno_t i=0; i < len; i++)
93  wgts[i] = val;
94  }
95  else if (how == roundRobin){
96  for (int i=0; i < 10; i++){
97  zscalar_t val = (i + 10)*scale;
98  for (zlno_t j=i; j < len; j += 10)
99  weights[j] = val;
100  }
101  }
102  else if (how == increasing){
103  zscalar_t val = scale + rank;
104  for (zlno_t i=0; i < len; i++)
105  wgts[i] = val;
106  }
107 
108  return weights;
109 }
110 
115 const RCP<tMVector_t> getMeshCoordinates(
116  const RCP<const Teuchos::Comm<int> > & comm,
117  zgno_t numGlobalCoords)
118 {
119  int rank = comm->getRank();
120  int nprocs = comm->getSize();
121 
122  double k = log(numGlobalCoords) / 3;
123  double xdimf = exp(k) + 0.5;
124  ssize_t xdim = static_cast<ssize_t>(floor(xdimf));
125  ssize_t ydim = xdim;
126  ssize_t zdim = numGlobalCoords / (xdim*ydim);
127  ssize_t num=xdim*ydim*zdim;
128  ssize_t diff = numGlobalCoords - num;
129  ssize_t newdiff = 0;
130 
131  while (diff > 0){
132  if (zdim > xdim && zdim > ydim){
133  zdim++;
134  newdiff = diff - (xdim*ydim);
135  if (newdiff < 0)
136  if (diff < -newdiff)
137  zdim--;
138  }
139  else if (ydim > xdim && ydim > zdim){
140  ydim++;
141  newdiff = diff - (xdim*zdim);
142  if (newdiff < 0)
143  if (diff < -newdiff)
144  ydim--;
145  }
146  else{
147  xdim++;
148  newdiff = diff - (ydim*zdim);
149  if (newdiff < 0)
150  if (diff < -newdiff)
151  xdim--;
152  }
153 
154  diff = newdiff;
155  }
156 
157  num=xdim*ydim*zdim;
158  diff = numGlobalCoords - num;
159  if (diff < 0)
160  diff /= -numGlobalCoords;
161  else
162  diff /= numGlobalCoords;
163 
164  if (rank == 0){
165  if (diff > .01)
166  std::cout << "Warning: Difference " << diff*100 << " percent" << std::endl;
167  std::cout << "Mesh size: " << xdim << "x" << ydim << "x" <<
168  zdim << ", " << num << " vertices." << std::endl;
169  }
170 
171  // Divide coordinates.
172 
173  ssize_t numLocalCoords = num / nprocs;
174  ssize_t leftOver = num % nprocs;
175  ssize_t gid0 = 0;
176 
177  if (rank <= leftOver)
178  gid0 = zgno_t(rank) * (numLocalCoords+1);
179  else
180  gid0 = (leftOver * (numLocalCoords+1)) +
181  ((zgno_t(rank) - leftOver) * numLocalCoords);
182 
183  if (rank < leftOver)
184  numLocalCoords++;
185 
186  ssize_t gid1 = gid0 + numLocalCoords;
187 
188  zgno_t *ids = new zgno_t [numLocalCoords];
189  if (!ids)
190  throw bad_alloc();
191  ArrayRCP<zgno_t> idArray(ids, 0, numLocalCoords, true);
192 
193  for (ssize_t i=gid0; i < gid1; i++)
194  *ids++ = zgno_t(i);
195 
196  RCP<const tMap_t> idMap = rcp(
197  new tMap_t(num, idArray.view(0, numLocalCoords), 0, comm));
198 
199  // Create a Tpetra::MultiVector of coordinates.
200 
201  zscalar_t *x = new zscalar_t [numLocalCoords*3];
202  if (!x)
203  throw bad_alloc();
204  ArrayRCP<zscalar_t> coordArray(x, 0, numLocalCoords*3, true);
205 
206  zscalar_t *y = x + numLocalCoords;
207  zscalar_t *z = y + numLocalCoords;
208 
209  zgno_t xStart = 0;
210  zgno_t yStart = 0;
211  zgno_t xyPlane = xdim*ydim;
212  zgno_t zStart = gid0 / xyPlane;
213  zgno_t rem = gid0 % xyPlane;
214  if (rem > 0){
215  yStart = rem / xdim;
216  xStart = rem % xdim;
217  }
218 
219  zlno_t next = 0;
220  for (zscalar_t zval=zStart; next < numLocalCoords && zval < zdim; zval++){
221  for (zscalar_t yval=yStart; next < numLocalCoords && yval < ydim; yval++){
222  for (zscalar_t xval=xStart; next < numLocalCoords && xval < xdim; xval++){
223  x[next] = xval;
224  y[next] = yval;
225  z[next] = zval;
226  next++;
227  }
228  xStart = 0;
229  }
230  yStart = 0;
231  }
232 
233  ArrayView<const zscalar_t> xArray(x, numLocalCoords);
234  ArrayView<const zscalar_t> yArray(y, numLocalCoords);
235  ArrayView<const zscalar_t> zArray(z, numLocalCoords);
236  ArrayRCP<ArrayView<const zscalar_t> > coordinates =
237  arcp(new ArrayView<const zscalar_t> [3], 0, 3);
238  coordinates[0] = xArray;
239  coordinates[1] = yArray;
240  coordinates[2] = zArray;
241 
242  ArrayRCP<const ArrayView<const zscalar_t> > constCoords =
243  coordinates.getConst();
244 
245  RCP<tMVector_t> meshCoords = rcp(new tMVector_t(
246  idMap, constCoords.view(0,3), 3));
247 
248  return meshCoords;
249 }
250 
251 
252 int main(int narg, char *arg[])
253 {
254  // MEMORY_CHECK(true, "Before initializing MPI");
255 
256  Tpetra::ScopeGuard tscope(&narg, &arg);
257  Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
258 
259  int rank = comm->getRank();
260  int nprocs = comm->getSize();
261 
262  MEMORY_CHECK(rank==0, "After initializing MPI");
263 
264  if (rank==0)
265  std::cout << "Number of processes: " << nprocs << std::endl;
266 
267  // Default values
268  double numGlobalCoords = 1000;
269  int numTestCuts = 1;
270  int nWeights = 0;
271  string timingType("no_timers");
272  string debugLevel("basic_status");
273  string memoryOn("memoryOn");
274  string memoryOff("memoryOff");
275  string memoryProcs("0");
276  bool doMemory=false;
277  int numGlobalParts = nprocs;
278 
279  CommandLineProcessor commandLine(false, true);
280  commandLine.setOption("size", &numGlobalCoords,
281  "Approximate number of global coordinates.");
282  commandLine.setOption("testCuts", &numTestCuts,
283  "Number of test cuts to make when looking for bisector.");
284  commandLine.setOption("numParts", &numGlobalParts,
285  "Number of parts (default is one per proc).");
286  commandLine.setOption("nWeights", &nWeights,
287  "Number of weights per coordinate, zero implies uniform weights.");
288 
289  string balanceCount("balance_object_count");
290  string balanceWeight("balance_object_weight");
291  string mcnorm1("multicriteria_minimize_total_weight");
292  string mcnorm2("multicriteria_balance_total_maximum");
293  string mcnorm3("multicriteria_minimize_maximum_weight");
294 
295  string objective(balanceWeight); // default
296 
297  string doc(balanceCount); doc.append(": ignore weights\n");
298 
299  doc.append(balanceWeight); doc.append(": balance on first weight\n");
300 
301  doc.append(mcnorm1);
302  doc.append(": given multiple weights, balance their total.\n");
303 
304  doc.append(mcnorm3);
305  doc.append(": given multiple weights, balance the maximum for each coordinate.\n");
306 
307  doc.append(mcnorm2);
308  doc.append(": given multiple weights, balance the L2 norm of the weights.\n");
309 
310  commandLine.setOption("objective", &objective, doc.c_str());
311 
312  commandLine.setOption("timers", &timingType,
313  "no_timers, micro_timers, macro_timers, both_timers, test_timers");
314 
315  commandLine.setOption("debug", &debugLevel,
316  "no_status, basic_status, detailed_status, verbose_detailed_status");
317 
318  commandLine.setOption(memoryOn.c_str(), memoryOff.c_str(), &doMemory,
319  "do memory profiling");
320 
321  commandLine.setOption("memoryProcs", &memoryProcs,
322  "list of processes that output memory usage");
323 
324  CommandLineProcessor::EParseCommandLineReturn rc =
325  commandLine.parse(narg, arg);
326 
327  if (rc != Teuchos::CommandLineProcessor::PARSE_SUCCESSFUL){
328  if (rc == Teuchos::CommandLineProcessor::PARSE_HELP_PRINTED){
329  if (rank==0)
330  std::cout << "PASS" << std::endl;
331  return 1;
332  }
333  else{
334  if (rank==0)
335  std::cout << "FAIL" << std::endl;
336  return 0;
337  }
338  }
339 
340  //MEMORY_CHECK(doMemory && rank==0, "After processing parameters");
341 
342  zgno_t globalSize = static_cast<zgno_t>(numGlobalCoords);
343 
344  RCP<tMVector_t> coordinates = getMeshCoordinates(comm, globalSize);
345  size_t numLocalCoords = coordinates->getLocalLength();
346 
347 #if 0
348  comm->barrier();
349  for (int p=0; p < nprocs; p++){
350  if (p==rank){
351  std::cout << "Rank " << rank << ", " << numLocalCoords << "coords" << std::endl;
352  const zscalar_t *x = coordinates->getData(0).getRawPtr();
353  const zscalar_t *y = coordinates->getData(1).getRawPtr();
354  const zscalar_t *z = coordinates->getData(2).getRawPtr();
355  for (zlno_t i=0; i < numLocalCoords; i++)
356  std::cout << " " << x[i] << " " << y[i] << " " << z[i] << std::endl;
357  }
358  std::cout.flush();
359  comm->barrier();
360  }
361 #endif
362 
363  Array<ArrayRCP<zscalar_t> > weights(nWeights);
364 
365  if (nWeights > 0){
366  int wt = 0;
367  zscalar_t scale = 1.0;
368  for (int i=0; i < nWeights; i++){
369  weights[i] =
370  makeWeights(comm, numLocalCoords, weightTypes(wt++), scale, rank);
371  if (wt == numWeightTypes){
372  wt = 0;
373  scale++;
374  }
375  }
376  }
377 
378  MEMORY_CHECK(doMemory && rank==0, "After creating input");
379 
380  // Create an input adapter.
381  const RCP<const tMap_t> &coordmap = coordinates->getMap();
382  ArrayView<const zgno_t> ids = coordmap->getLocalElementList();
383  const zgno_t *globalIds = ids.getRawPtr();
384 
385  size_t localCount = coordinates->getLocalLength();
386  typedef Zoltan2::BasicVectorAdapter<tMVector_t> inputAdapter_t;
387  RCP<inputAdapter_t> ia;
388 
389  if (nWeights == 0){
390  ia = rcp(new inputAdapter_t (localCount, globalIds,
391  coordinates->getData(0).getRawPtr(), coordinates->getData(1).getRawPtr(),
392  coordinates->getData(2).getRawPtr(), 1,1,1));
393  }
394  else{
395  vector<const zscalar_t *> values(3);
396  for (int i=0; i < 3; i++)
397  values[i] = coordinates->getData(i).getRawPtr();
398  vector<int> valueStrides(0); // implies stride is one
399  vector<const zscalar_t *> weightPtrs(nWeights);
400  for (int i=0; i < nWeights; i++)
401  weightPtrs[i] = weights[i].getRawPtr();
402  vector<int> weightStrides(0); // implies stride is one
403 
404  ia = rcp(new inputAdapter_t (localCount, globalIds,
405  values, valueStrides, weightPtrs, weightStrides));
406  }
407 
408  MEMORY_CHECK(doMemory && rank==0, "After creating input adapter");
409 
410  // Parameters
411 
412  Teuchos::ParameterList params;
413 
414  if (timingType != "no_timers"){
415  params.set("timer_output_stream" , "std::cout");
416  params.set("timer_type" , timingType);
417  }
418 
419  if (doMemory){
420  params.set("memory_output_stream" , "std::cout");
421  params.set("memory_procs" , memoryProcs);
422  }
423 
424  params.set("debug_output_stream" , "std::cerr");
425  params.set("debug_procs" , "0");
426 
427  if (debugLevel != "basic_status"){
428  params.set("debug_level" , debugLevel);
429  }
430 
431  params.set("algorithm", "rcb");
432  params.set("partitioning_objective", objective);
433  double tolerance = 1.1;
434  params.set("imbalance_tolerance", tolerance );
435 
436  if (numGlobalParts != nprocs)
437  params.set("num_global_parts" , numGlobalParts);
438 
439  if (rank==0){
440  std::cout << "Number of parts: " << numGlobalParts << std::endl;
441  }
442 
443  // Create a problem, solve it, and display the quality.
444 
445  Zoltan2::PartitioningProblem<inputAdapter_t> problem(&(*ia), &params);
446 
447  problem.solve();
448 
449  comm->barrier();
450 
451  problem.printTimers();
452 
453  comm->barrier();
454 
455  if (rank == 0){
456  std::cout << "PASS" << std::endl;
457  }
458 
459  return 0;
460 }
461 
const RCP< tMVector_t > getMeshCoordinates(const RCP< const Teuchos::Comm< int > > &comm, zgno_t numGlobalCoords)
Create a mesh of approximately the desired size.
static ArrayRCP< ArrayRCP< zscalar_t > > weights
int main(int narg, char **arg)
Definition: coloring1.cpp:199
Defines the PartitioningSolution class.
common code used by tests
ArrayRCP< zscalar_t > makeWeights(const RCP< const Teuchos::Comm< int > > &comm, zlno_t len, weightTypes how, zscalar_t scale, int rank)
static RCP< tMVector_t > coordinates
BasicVectorAdapter represents a vector (plus optional weights) supplied by the user as pointers to st...
void printTimers() const
Return the communicator passed to the problem.
Tpetra::Map::local_ordinal_type zlno_t
Tpetra::Map< zlno_t, zgno_t, znode_t > tMap_t
PartitioningProblem sets up partitioning problems for the user.
Defines the PartitioningProblem class.
float zscalar_t
Defines the BasicVectorAdapter class.
Tpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > tMVector_t
Tpetra::Map::global_ordinal_type zgno_t
void solve(bool updateInputData=true)
Direct the problem to create a solution.
#define MEMORY_CHECK(iPrint, msg)