Zoltan2
Loading...
Searching...
No Matches
blockTest.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
49
50#include <Teuchos_DefaultComm.hpp>
51#include <Teuchos_CommandLineProcessor.hpp>
52#include <Teuchos_ParameterList.hpp>
53
54
55int main(int narg, char **arg)
56{
57 Tpetra::ScopeGuard tscope(&narg, &arg);
58 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
59
60 int fail=0, gfail=0;
61
62 int rank = comm->getRank();
63 int nprocs = comm->getSize();
64
65 // Construct the user input
66 int numGlobalIdentifiers = 100;
67 int numMyIdentifiers = numGlobalIdentifiers / nprocs;
68 if (rank < numGlobalIdentifiers % nprocs)
69 numMyIdentifiers += 1;
70
71 zgno_t myBaseId = zgno_t(numGlobalIdentifiers * rank);
72
73 zgno_t *myIds = new zgno_t[numMyIdentifiers];
74 zscalar_t *myWeights = new zscalar_t[numMyIdentifiers];
75
76 if (!myIds || !myWeights){
77 fail = 1;
78 }
79
80 gfail = globalFail(*comm, fail);
81
82 if (gfail){
83 if (rank==0){
84 std::cout << "Memory allocation failure" << std::endl;
85 std::cout << "FAIL" << std::endl;
86 }
87 return 1;
88 }
89
90 zscalar_t origsumwgts = 0;
91 for (int i=0; i < numMyIdentifiers; i++){
92 myIds[i] = myBaseId+i;
93 myWeights[i] = rank%3 + 1;
94 origsumwgts += myWeights[i];
95 }
96
97 // Some output
98 int *origcnt = new int[nprocs];
99 zscalar_t *origwgts = new zscalar_t[nprocs];
100 Teuchos::gather<int, int>(&numMyIdentifiers, 1, origcnt, 1, 0, *comm);
101 Teuchos::gather<int, zscalar_t>(&origsumwgts, 1, origwgts, 1, 0, *comm);
102 if (rank == 0) {
103 std::cout << "BEFORE PART CNTS: ";
104 for (int i = 0; i < nprocs; i++)
105 std::cout << origcnt[i] << " ";
106 std::cout << std::endl;
107 std::cout << "BEFORE PART WGTS: ";
108 for (int i = 0; i < nprocs; i++)
109 std::cout << origwgts[i] << " ";
110 std::cout << std::endl;
111 }
112 delete [] origcnt;
113 delete [] origwgts;
114
115 // Building Zoltan2 adapters
116 std::vector<const zscalar_t *> weightValues;
117 std::vector<int> weightStrides; // default is one
118 weightValues.push_back(const_cast<const zscalar_t *>(myWeights));
119
123 typedef adapter_t::part_t part_t;
124
125 adapter_t *adapter =
126 new adapter_t(zlno_t(numMyIdentifiers),myIds,weightValues,weightStrides);
127
128 // Set up the parameters and problem
129 bool useWeights = true;
130 Teuchos::CommandLineProcessor cmdp (false, false);
131 cmdp.setOption("weights", "no-weights", &useWeights,
132 "Indicated whether to use identifier weights in partitioning");
133 cmdp.parse(narg, arg);
134
135 Teuchos::ParameterList params("test parameters");
136 //params.set("compute_metrics", true); // bool parameter
137 params.set("num_global_parts", nprocs);
138 params.set("algorithm", "block");
139 params.set("partitioning_approach", "partition");
140 if (!useWeights) params.set("partitioning_objective", "balance_object_count");
141
142 Zoltan2::PartitioningProblem<adapter_t> problem(adapter, &params);
143
144 problem.solve();
145
147
148 // create metric object
149
150 quality_t *metricObject = new quality_t(adapter, &params, comm, &solution);
151
152 // Some output
153 zscalar_t *totalWeight = new zscalar_t [nprocs];
154 zscalar_t *sumWeight = new zscalar_t [nprocs];
155 memset(totalWeight, 0, nprocs * sizeof(zscalar_t));
156 int *totalCnt = new int [nprocs];
157 int *sumCnt = new int [nprocs];
158 memset(totalCnt, 0, nprocs * sizeof(int));
159
160 const part_t *partList = solution.getPartListView();
161 zscalar_t libImbalance = metricObject->getWeightImbalance(0);
162 delete metricObject;
163 delete adapter;
164
165 for (int i=0; i < numMyIdentifiers; i++){
166 totalCnt[partList[i]]++;
167 totalWeight[partList[i]] += myWeights[i];
168 }
169
170 Teuchos::reduceAll<int, int>(*comm, Teuchos::REDUCE_SUM, nprocs,
171 totalCnt, sumCnt);
172 Teuchos::reduceAll<int, zscalar_t>(*comm, Teuchos::REDUCE_SUM, nprocs,
173 totalWeight, sumWeight);
174
175 double epsilon = 10e-6;
176
177 if (rank == 0){
178 std::cout << "AFTER PART CNTS: ";
179 for (int i=0; i < nprocs; i++)
180 std::cout << sumCnt[i] << " ";
181 std::cout << std::endl;
182
183 zscalar_t total = 0;
184 std::cout << "AFTER PART WGTS: ";
185 for (int i=0; i < nprocs; i++){
186 std::cout << sumWeight[i] << " ";
187 total += sumWeight[i];
188 }
189 std::cout << std::endl;
190
191 zscalar_t avg = total / zscalar_t(nprocs);
192
193 zscalar_t imbalance = -1.0;
194
195 for (int i=0; i < nprocs; i++){
196 zscalar_t imb = 0;
197 if (sumWeight[i] > avg)
198 imb = (sumWeight[i] - avg) / avg;
199 else
200 imb = (avg - sumWeight[i]) / avg;
201
202 if (imb > imbalance)
203 imbalance = imb;
204 }
205 imbalance += 1.0;
206
207 std::cout << "Computed imbalance: " << imbalance << std::endl;
208 std::cout << "Library's imbalance: " << libImbalance << std::endl;
209
210 double err;
211 if (imbalance > libImbalance)
212 err = imbalance - libImbalance;
213 else
214 err = libImbalance - imbalance;
215
216 if (err > epsilon)
217 fail = 1;
218 }
219 else{
220 fail = 0;
221 }
222
223 gfail = globalFail(*comm, fail);
224
225 if (gfail){
226 if (rank==0){
227 std::cout << "failure in solution's imbalance data" << std::endl;
228 std::cout << "FAIL" << std::endl;
229 }
230 return 1;
231 }
232
233 if (rank==0)
234 std::cout << "PASS" << std::endl;
235
236 delete [] myWeights;
237 delete [] myIds;
238 delete [] sumCnt;
239 delete [] totalCnt;
240 delete [] sumWeight;
241 delete [] totalWeight;
242
243 return 0;
244}
int globalFail(const Comm< int > &comm, int fail)
Defines the BasicIdentifierAdapter class.
Defines the PartitioningProblem class.
Defines the PartitioningSolution class.
common code used by tests
float zscalar_t
Tpetra::Map ::local_ordinal_type zlno_t
Tpetra::Map ::global_ordinal_type zgno_t
int main()
This class represents a collection of global Identifiers and their associated weights,...
A simple class that can be the User template argument for an InputAdapter.
A class that computes and returns quality metrics.
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.
A PartitioningSolution is a solution to a partitioning problem.
static const std::string fail
#define epsilon
Definition nd.cpp:82
SparseMatrixAdapter_t::part_t part_t
Zoltan2::EvaluatePartition< matrixAdapter_t > quality_t