Zoltan2
Loading...
Searching...
No Matches
fix4785.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
56#include <Tpetra_Map.hpp>
57#include <vector>
58#include <cstdlib>
59
60typedef Tpetra::Map<> myMap_t;
61typedef myMap_t::local_ordinal_type myLocalId_t;
62typedef myMap_t::global_ordinal_type myGlobalId_t;
63typedef double myScalar_t;
64
66template <typename A, typename B>
67void copyArrays(size_t n, A *&a, B *b)
68{
69 a = new A[n];
70 for (size_t i = 0; i < n; i++) a[i] = b[i];
71}
72
73template <typename A, typename B>
74void copyArrays(size_t n, A *&a, A *b)
75{
76 a = b;
77}
78
80template <typename A, typename B>
81void freeArrays(A *&a, B *b)
82{
83 delete [] a;
84}
85
86template <typename A, typename B>
87void freeArrays(A *&a, A *b)
88{
89 // no delete needed since only copied pointer
90}
91
94 const Teuchos::RCP<const Teuchos::Comm<int> > &comm,
95 Teuchos::ParameterList &params,
96 size_t localCount,
97 myGlobalId_t *globalIds,
98 myScalar_t *coords,
99 int &nFail
100)
101{
102 typedef Tpetra::Map<>::node_type myNode_t;
103
105 myNode_t> myTypes;
106 typedef Zoltan2::BasicVectorAdapter<myTypes> inputAdapter_t;
108
109 myScalar_t *sx = coords;
110 myScalar_t *sy = sx + localCount;
111 myScalar_t *sz = sy + localCount;
112
113 inputAdapter_t *ia = new inputAdapter_t(localCount,globalIds,sx,sy,sz,1,1,1);
114
117
118 problem->solve();
119
120 quality_t *metricObject = new quality_t(ia, &params, comm,
121 &problem->getSolution());
122 if (comm->getRank() == 0){
123
124 metricObject->printMetrics(std::cout);
125
126 double imb = metricObject->getObjectCountImbalance();
127 if (imb <= 1.01) // Should get perfect balance
128 std::cout << "no weights -- balance satisfied: " << imb << std::endl;
129 else {
130 std::cout << "no weights -- balance failure: " << imb << std::endl;
131 nFail++;
132 }
133 std::cout << std::endl;
134 }
135
136 delete metricObject;
137 delete problem;
138 delete ia;
139}
140
143 const Teuchos::RCP<const Teuchos::Comm<int> > &comm,
144 Teuchos::ParameterList &params,
145 size_t localCount,
146 myGlobalId_t *globalIds,
147 myScalar_t *coords,
149 int &nFail
150)
151{
152 typedef Tpetra::Map<>::node_type myNode_t;
153
155 myNode_t> myTypes;
156 typedef Zoltan2::BasicVectorAdapter<myTypes> inputAdapter_t;
158
159 std::vector<const myScalar_t *> coordVec(3);
160 std::vector<int> coordStrides(3);
161
162 coordVec[0] = coords; coordStrides[0] = 1;
163 coordVec[1] = coords + localCount; coordStrides[1] = 1;
164 coordVec[2] = coords + localCount + localCount; coordStrides[2] = 1;
165
166 std::vector<const myScalar_t *> weightVec(1);
167 std::vector<int> weightStrides(1);
168
169 weightVec[0] = weights; weightStrides[0] = 1;
170
171 inputAdapter_t *ia=new inputAdapter_t(localCount, globalIds, coordVec,
172 coordStrides,weightVec,weightStrides);
173
176
177 problem->solve();
178
179 quality_t *metricObject = new quality_t(ia, &params, comm,
180 &problem->getSolution());
181 if (comm->getRank() == 0){
182
183 metricObject->printMetrics(std::cout);
184
185 double imb = metricObject->getWeightImbalance(0);
186 if (imb <= 1.01)
187 std::cout << "weighted -- balance satisfied " << imb << std::endl;
188 else {
189 std::cout << "weighted -- balance failed " << imb << std::endl;
190 nFail++;
191 }
192 std::cout << std::endl;
193 }
194
195 delete metricObject;
196 delete problem;
197 delete ia;
198}
199
201int main(int narg, char *arg[])
202{
203 Tpetra::ScopeGuard scope(&narg, &arg);
204 const Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
205 int me = comm->getRank();
206 int np = comm->getSize();
207 int nFail = 0;
208
210 // Create parameters for an MJ problem
211
212 Teuchos::ParameterList params("test params");
213 params.set("debug_level", "basic_status");
214 params.set("error_check_level", "debug_mode_assertions");
215
216 params.set("algorithm", "multijagged");
217 params.set("num_global_parts", 64);
218
220 // Create input data.
221
222 const size_t N = 9000000;
223 size_t maxLocalCount = N / np; // biggest test we'll run
224
225 // Create coordinates that range from 0 to 999
226 const int dim = 3;
227 myScalar_t *coords = new myScalar_t[dim * maxLocalCount];
228
229 srand(me);
230 for (size_t i=0; i < maxLocalCount*dim; i++)
231 coords[i] = myScalar_t(0);
232
233 // Create weights for the coordinates
234 myScalar_t *weights = new myScalar_t[maxLocalCount];
235 for (size_t i=0; i < maxLocalCount; i++)
236 weights[i] = myScalar_t(1+me);
237
238 // Allocate space for global IDs; they will be generated below
239 myGlobalId_t *globalIds = new myGlobalId_t[maxLocalCount];
240
241 size_t localCount = N / np;
242
243 // Create consecutive global ids for the coordinates
244 myGlobalId_t offset = me * localCount;
245 for (size_t i=0; i < localCount; i++)
246 globalIds[i] = offset++;
247
248 if (me == 0) {
249 std::cout << "---------------------------------------------" << std::endl;
250 std::cout << "myGlobalId_t = " << typeid(offset).name()
251 << " " << sizeof(offset)
252 << "; localCount = " << localCount
253 << "; globalCount = " << np * localCount << std::endl;
254 }
255
257 // Test one: No weights
258 if (me == 0) std::cout << "Test: no weights, scalar = double" << std::endl;
259 test_no_weights(comm, params, localCount, globalIds, coords, nFail);
260
262 // Test two: weighted
263 if (me == 0) std::cout << "Test: weights, scalar = double" << std::endl;
264 test_weights(comm, params, localCount, globalIds, coords, weights, nFail);
265
266 // Early exit when failure is detected
267 int gnFail;
268 Teuchos::reduceAll<int, int>(*comm, Teuchos::REDUCE_MAX, 1, &nFail, &gnFail);
269
270 delete [] weights;
271 delete [] coords;
272 delete [] globalIds;
273
274 if (me == 0) {
275 if (nFail == 0) std::cout << "PASS" << std::endl;
276 else std::cout << "FAIL: " << nFail << " tests failed" << std::endl;
277 }
278
279 return 0;
280}
281
Defines the BasicVectorAdapter class.
Traits for application input objects.
Defines the PartitioningProblem class.
Defines the PartitioningSolution class.
int main()
A simple class that can be the User template argument for an InputAdapter.
BasicVectorAdapter represents a vector (plus optional weights) supplied by the user as pointers to st...
A class that computes and returns quality metrics.
void printMetrics(std::ostream &os) const
Print all metrics.
scalar_t getWeightImbalance(int weightIndex) const
Return the imbalance for the requested weight.
scalar_t getObjectCountImbalance() const
Return the object count imbalance.
PartitioningProblem sets up partitioning problems for the user.
const PartitioningSolution< Adapter > & getSolution()
Get the solution to the problem.
void solve(bool updateInputData=true)
Direct the problem to create a solution.
myMap_t::local_ordinal_type myLocalId_t
Definition: fix4785.cpp:61
double myScalar_t
Definition: fix4785.cpp:63
void freeArrays(A *&a, B *b)
Definition: fix4785.cpp:81
void test_no_weights(const Teuchos::RCP< const Teuchos::Comm< int > > &comm, Teuchos::ParameterList &params, size_t localCount, myGlobalId_t *globalIds, myScalar_t *coords, int &nFail)
Definition: fix4785.cpp:93
Tpetra::Map myMap_t
Definition: fix4785.cpp:60
void copyArrays(size_t n, A *&a, B *b)
Definition: fix4785.cpp:67
void test_weights(const Teuchos::RCP< const Teuchos::Comm< int > > &comm, Teuchos::ParameterList &params, size_t localCount, myGlobalId_t *globalIds, myScalar_t *coords, myScalar_t *weights, int &nFail)
Definition: fix4785.cpp:142
myMap_t::global_ordinal_type myGlobalId_t
Definition: fix4785.cpp:62
static ArrayRCP< ArrayRCP< zscalar_t > > weights
Zoltan2::EvaluatePartition< matrixAdapter_t > quality_t