Zoltan2
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
mj_int_coordinates.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 
54 #include <Zoltan2_InputTraits.hpp>
55 #include <Tpetra_Map.hpp>
56 #include <vector>
57 #include <cstdlib>
58 
59 int main(int narg, char *arg[])
60 {
61  Tpetra::ScopeGuard scope(&narg, &arg);
62  const Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
63  int rank = comm->getRank();
64  int nprocs = comm->getSize();
65  int nFail = 0;
66 
67  typedef Tpetra::Map<> Map_t;
68  typedef Map_t::local_ordinal_type localId_t;
69  typedef Map_t::global_ordinal_type globalId_t;
70  typedef int int_scalar_t; // This is the case we are testing here.
71 
73  typedef Zoltan2::BasicVectorAdapter<myTypes> inputAdapter_t;
74 
76 
78  // Create input data.
79 
80  size_t localCount = 40;
81  int dim = 3;
82 
83  // Create coordinates that range from 0 to 9
84  int_scalar_t *coords = new int_scalar_t [dim * localCount];
85  int_scalar_t *x = coords;
86  int_scalar_t *y = x + localCount;
87  int_scalar_t *z = y + localCount;
88 
89  srand(rank);
90  for (size_t i=0; i < localCount*dim; i++)
91  coords[i] = int_scalar_t(rand() % 10);
92 
93  // Create global ids for the coordinates.
94  globalId_t *globalIds = new globalId_t [localCount];
95  globalId_t offset = rank * localCount;
96  for (size_t i=0; i < localCount; i++) globalIds[i] = offset++;
97 
99  // Create parameters for an MJ problem
100 
101  Teuchos::ParameterList params("test params");
102  params.set("debug_level", "basic_status");
103  params.set("error_check_level", "debug_mode_assertions");
104 
105  params.set("algorithm", "multijagged");
106  params.set("num_global_parts", nprocs+1);
107 
109  // Test one: No weights
110 
111  inputAdapter_t *ia1 = new inputAdapter_t(localCount,globalIds,x,y,z,1,1,1);
112 
115 
116  problem1->solve();
117 
118  quality_t *metricObject1 = new quality_t(ia1, &params, comm,
119  &problem1->getSolution());
120  if (rank == 0){
121 
122  metricObject1->printMetrics(std::cout);
123 
124  double imb = metricObject1->getObjectCountImbalance();
125  if (imb <= 1.01) // Should get perfect balance
126  std::cout << "no weights -- balance satisfied: " << imb << std::endl;
127  else {
128  std::cout << "no weights -- balance failure: " << imb << std::endl;
129  nFail++;
130  }
131  std::cout << std::endl;
132  }
133  delete metricObject1;
134  delete problem1;
135  delete ia1;
136 
138  // Test two: weighted
139  // Create a Zoltan2 input adapter that includes weights.
140 
141  int_scalar_t *weights = new int_scalar_t [localCount];
142  for (size_t i=0; i < localCount; i++) weights[i] = 1 + int_scalar_t(rank);
143 
144  std::vector<const int_scalar_t *>coordVec(3);
145  std::vector<int> coordStrides(3);
146 
147  coordVec[0] = x; coordStrides[0] = 1;
148  coordVec[1] = y; coordStrides[1] = 1;
149  coordVec[2] = z; coordStrides[2] = 1;
150 
151  std::vector<const int_scalar_t *>weightVec(1);
152  std::vector<int> weightStrides(1);
153 
154  weightVec[0] = weights; weightStrides[0] = 1;
155 
156  inputAdapter_t *ia2=new inputAdapter_t(localCount, globalIds, coordVec,
157  coordStrides,weightVec,weightStrides);
158 
161 
162  problem2->solve();
163 
164  quality_t *metricObject2 = new quality_t(ia2, &params, comm,
165  &problem2->getSolution());
166  if (rank == 0){
167 
168  metricObject2->printMetrics(std::cout);
169 
170  double imb = metricObject2->getWeightImbalance(0);
171  if (imb <= 1.01)
172  std::cout << "weighted -- balance satisfied " << imb << std::endl;
173  else {
174  std::cout << "weighted -- balance failed " << imb << std::endl;
175  nFail++;
176  }
177  std::cout << std::endl;
178  }
179  delete metricObject2;
180 
181  if (weights) delete [] weights;
182  if (coords) delete [] coords;
183  if (globalIds) delete [] globalIds;
184  delete problem2;
185  delete ia2;
186 
187  if (rank == 0) {
188  if (nFail == 0) std::cout << "PASS" << std::endl;
189  else std::cout << "FAIL: " << nFail << " tests failed" << std::endl;
190  }
191 
192  return 0;
193 }
194 
int main(int narg, char *arg[])
A simple class that can be the User template argument for an InputAdapter.
static ArrayRCP< ArrayRCP< zscalar_t > > weights
Defines the PartitioningSolution class.
BasicVectorAdapter represents a vector (plus optional weights) supplied by the user as pointers to st...
Traits for application input objects.
const PartitioningSolution< Adapter > & getSolution()
Get the solution to the problem.
PartitioningProblem sets up partitioning problems for the user.
Defines the PartitioningProblem class.
A class that computes and returns quality metrics.
Defines the BasicVectorAdapter class.
void solve(bool updateInputData=true)
Direct the problem to create a solution.