Zoltan2
 All Namespaces Files Functions Variables Pages
rcb_C.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 
50 #include <Zoltan2_PartitioningSolution.hpp>
51 #include <Zoltan2_PartitioningProblem.hpp>
52 #include <Zoltan2_BasicVectorAdapter.hpp>
53 #include <Zoltan2_InputTraits.hpp>
54 #include <Tpetra_Map.hpp>
55 #include <vector>
56 #include <cstdlib>
57 
62 int main(int argc, char *argv[])
63 {
64  Tpetra::ScopeGuard tscope(&argc, &argv);
65  Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
66 
67  int rank = comm->getRank();
68  int nprocs = comm->getSize();
69 
70  // For convenience, we'll use the Tpetra defaults for local/global ID types
71  // Users can substitute their preferred local/global ID types
72  typedef Tpetra::Map<> Map_t;
73  typedef Map_t::local_ordinal_type localId_t;
74  typedef Map_t::global_ordinal_type globalId_t;
75 
76  typedef Tpetra::Details::DefaultTypes::scalar_type scalar_t;
77  typedef Zoltan2::BasicUserTypes<scalar_t, localId_t, globalId_t> myTypes;
78 
79  // TODO explain
80  typedef Zoltan2::BasicVectorAdapter<myTypes> inputAdapter_t;
81  typedef Zoltan2::EvaluatePartition<inputAdapter_t> quality_t;
82  typedef inputAdapter_t::part_t part_t;
83 
85  // Create input data.
86 
87  size_t localCount = 40;
88  int dim = 3;
89 
90  scalar_t *coords = new scalar_t [dim * localCount];
91 
92  scalar_t *x = coords;
93  scalar_t *y = x + localCount;
94  scalar_t *z = y + localCount;
95 
96  // Create coordinates that range from 0 to 10.0
97 
98  srand(rank);
99  scalar_t scalingFactor = 10.0 / RAND_MAX;
100 
101  for (size_t i=0; i < localCount*dim; i++){
102  coords[i] = scalar_t(rand()) * scalingFactor;
103  }
104 
105  // Create global ids for the coordinates.
106 
107  globalId_t *globalIds = new globalId_t [localCount];
108  globalId_t offset = rank * localCount;
109 
110  for (size_t i=0; i < localCount; i++)
111  globalIds[i] = offset++;
112 
114  // Create parameters for an RCB problem
115 
116  double tolerance = 1.1;
117 
118  if (rank == 0)
119  std::cout << "Imbalance tolerance is " << tolerance << std::endl;
120 
121  Teuchos::ParameterList params("test params");
122  params.set("debug_level", "basic_status");
123  params.set("debug_procs", "0");
124  params.set("error_check_level", "debug_mode_assertions");
125 
126  params.set("algorithm", "rcb");
127  params.set("imbalance_tolerance", tolerance );
128  params.set("num_global_parts", nprocs);
129 
132  // A simple problem with no weights.
135 
136  // Create a Zoltan2 input adapter for this geometry. TODO explain
137 
138  inputAdapter_t *ia1 = new inputAdapter_t(localCount,globalIds,x,y,z,1,1,1);
139 
140  // Create a Zoltan2 partitioning problem
141 
142  Zoltan2::PartitioningProblem<inputAdapter_t> *problem1 =
143  new Zoltan2::PartitioningProblem<inputAdapter_t>(ia1, &params);
144 
145  // Solve the problem
146 
147  problem1->solve();
148 
149  // create metric object where communicator is Teuchos default
150 
151  quality_t *metricObject1 = new quality_t(ia1, &params, //problem1->getComm(),
152  &problem1->getSolution());
153  // Check the solution.
154 
155  if (rank == 0) {
156  metricObject1->printMetrics(std::cout);
157  }
158 
159  if (rank == 0){
160  scalar_t imb = metricObject1->getObjectCountImbalance();
161  if (imb <= tolerance)
162  std::cout << "pass: " << imb << std::endl;
163  else
164  std::cout << "fail: " << imb << std::endl;
165  std::cout << std::endl;
166  }
167  delete metricObject1;
168 
171  // Try a problem with weights
174 
175  scalar_t *weights = new scalar_t [localCount];
176  for (size_t i=0; i < localCount; i++){
177  weights[i] = 1.0 + scalar_t(rank) / scalar_t(nprocs);
178  }
179 
180  // Create a Zoltan2 input adapter that includes weights.
181 
182  std::vector<const scalar_t *>coordVec(2);
183  std::vector<int> coordStrides(2);
184 
185  coordVec[0] = x; coordStrides[0] = 1;
186  coordVec[1] = y; coordStrides[1] = 1;
187 
188  std::vector<const scalar_t *>weightVec(1);
189  std::vector<int> weightStrides(1);
190 
191  weightVec[0] = weights; weightStrides[0] = 1;
192 
193  inputAdapter_t *ia2=new inputAdapter_t(localCount, globalIds, coordVec,
194  coordStrides,weightVec,weightStrides);
195 
196  // Create a Zoltan2 partitioning problem
197 
198  Zoltan2::PartitioningProblem<inputAdapter_t> *problem2 =
199  new Zoltan2::PartitioningProblem<inputAdapter_t>(ia2, &params);
200 
201  // Solve the problem
202 
203  problem2->solve();
204 
205  // create metric object for MPI builds
206 
207 #ifdef HAVE_ZOLTAN2_MPI
208  quality_t *metricObject2 = new quality_t(ia2, &params, //problem2->getComm()
209  MPI_COMM_WORLD,
210  &problem2->getSolution());
211 #else
212  quality_t *metricObject2 = new quality_t(ia2, &params, problem2->getComm(),
213  &problem2->getSolution());
214 #endif
215  // Check the solution.
216 
217  if (rank == 0) {
218  metricObject2->printMetrics(std::cout);
219  }
220 
221  if (rank == 0){
222  scalar_t imb = metricObject2->getWeightImbalance(0);
223  if (imb <= tolerance)
224  std::cout << "pass: " << imb << std::endl;
225  else
226  std::cout << "fail: " << imb << std::endl;
227  std::cout << std::endl;
228  }
229  delete metricObject2;
230 
231  if (localCount > 0){
232  delete [] weights;
233  weights = NULL;
234  }
235 
238  // Try a problem with multiple weights.
241 
242  // Add to the parameters the multicriteria objective.
243 
244  params.set("partitioning_objective", "multicriteria_minimize_total_weight");
245 
246  // Create the new weights.
247 
248  weights = new scalar_t [localCount*3];
249  srand(rank);
250 
251  for (size_t i=0; i < localCount*3; i+=3){
252  weights[i] = 1.0 + rank / nprocs; // weight idx 1
253  weights[i+1] = rank<nprocs/2 ? 1 : 2; // weight idx 2
254  weights[i+2] = rand()/RAND_MAX +.5; // weight idx 3
255  }
256 
257  // Create a Zoltan2 input adapter with these weights.
258 
259  weightVec.resize(3);
260  weightStrides.resize(3);
261 
262  weightVec[0] = weights; weightStrides[0] = 3;
263  weightVec[1] = weights+1; weightStrides[1] = 3;
264  weightVec[2] = weights+2; weightStrides[2] = 3;
265 
266  inputAdapter_t *ia3=new inputAdapter_t(localCount, globalIds, coordVec,
267  coordStrides,weightVec,weightStrides);
268 
269  // Create a Zoltan2 partitioning problem.
270 
271  Zoltan2::PartitioningProblem<inputAdapter_t> *problem3 =
272  new Zoltan2::PartitioningProblem<inputAdapter_t>(ia3, &params);
273 
274  // Solve the problem
275 
276  problem3->solve();
277 
278  // create metric object where Teuchos communicator is specified
279 
280  quality_t *metricObject3 = new quality_t(ia3, &params, problem3->getComm(),
281  &problem3->getSolution());
282  // Check the solution.
283 
284  if (rank == 0) {
285  metricObject3->printMetrics(std::cout);
286  }
287 
288  if (rank == 0){
289  scalar_t imb = metricObject3->getWeightImbalance(0);
290  if (imb <= tolerance)
291  std::cout << "pass: " << imb << std::endl;
292  else
293  std::cout << "fail: " << imb << std::endl;
294  std::cout << std::endl;
295  }
296  delete metricObject3;
297 
299  // Try the other multicriteria objectives.
300 
301  bool dataHasChanged = false; // default is true
302 
303  params.set("partitioning_objective", "multicriteria_minimize_maximum_weight");
304  problem3->resetParameters(&params);
305  problem3->solve(dataHasChanged);
306 
307  // Solution changed!
308 
309  metricObject3 = new quality_t(ia3, &params, problem3->getComm(),
310  &problem3->getSolution());
311  if (rank == 0){
312  metricObject3->printMetrics(std::cout);
313  scalar_t imb = metricObject3->getWeightImbalance(0);
314  if (imb <= tolerance)
315  std::cout << "pass: " << imb << std::endl;
316  else
317  std::cout << "fail: " << imb << std::endl;
318  std::cout << std::endl;
319  }
320  delete metricObject3;
321 
322  params.set("partitioning_objective", "multicriteria_balance_total_maximum");
323  problem3->resetParameters(&params);
324  problem3->solve(dataHasChanged);
325 
326  // Solution changed!
327 
328  metricObject3 = new quality_t(ia3, &params, problem3->getComm(),
329  &problem3->getSolution());
330  if (rank == 0){
331  metricObject3->printMetrics(std::cout);
332  scalar_t imb = metricObject3->getWeightImbalance(0);
333  if (imb <= tolerance)
334  std::cout << "pass: " << imb << std::endl;
335  else
336  std::cout << "fail: " << imb << std::endl;
337  std::cout << std::endl;
338  }
339  delete metricObject3;
340 
341  if (localCount > 0){
342  delete [] weights;
343  weights = NULL;
344  }
345 
348  // Using part sizes, ask for some parts to be empty.
351 
352  // Change the number of parts to twice the number of processes to
353  // ensure that we have more than one global part.
354 
355  params.set("num_global_parts", nprocs*2);
356 
357  // Using the initial problem that did not have any weights, reset
358  // parameter list, and give it some part sizes.
359 
360  problem1->resetParameters(&params);
361 
362  part_t partIds[2];
363  scalar_t partSizes[2];
364 
365  partIds[0] = rank*2; partSizes[0] = 0;
366  partIds[1] = rank*2+1; partSizes[1] = 1;
367 
368  problem1->setPartSizes(2, partIds, partSizes);
369 
370  // Solve the problem. The argument "dataHasChanged" indicates
371  // that we have not changed the input data, which allows the problem
372  // so skip some work when re-solving.
373 
374  dataHasChanged = false;
375 
376  problem1->solve(dataHasChanged);
377 
378  // Obtain the solution
379 
380  const Zoltan2::PartitioningSolution<inputAdapter_t> &solution4 =
381  problem1->getSolution();
382 
383  // Check it. Part sizes should all be odd.
384 
385  const part_t *partAssignments = solution4.getPartListView();
386 
387  int numInEmptyParts = 0;
388  for (size_t i=0; i < localCount; i++){
389  if (partAssignments[i] % 2 == 0)
390  numInEmptyParts++;
391  }
392 
393  if (rank == 0)
394  std::cout << "Request that " << nprocs << " parts be empty." <<std::endl;
395 
396  // Solution changed!
397 
398  metricObject1 = new quality_t(ia1, &params, //problem1->getComm(),
399  &problem1->getSolution());
400  // Check the solution.
401 
402  if (rank == 0) {
403  metricObject1->printMetrics(std::cout);
404  }
405 
406  if (rank == 0){
407  scalar_t imb = metricObject1->getObjectCountImbalance();
408  if (imb <= tolerance)
409  std::cout << "pass: " << imb << std::endl;
410  else
411  std::cout << "fail: " << imb << std::endl;
412  std::cout << std::endl;
413  }
414  delete metricObject1;
415 
416  if (coords)
417  delete [] coords;
418 
419  if (globalIds)
420  delete [] globalIds;
421 
422  delete problem1;
423  delete ia1;
424  delete problem2;
425  delete ia2;
426  delete problem3;
427  delete ia3;
428 
429  if (rank == 0)
430  std::cout << "PASS" << std::endl;
431 }
432 
int main(int argc, char *argv[])
Definition: block.cpp:60