Zoltan2
Loading...
Searching...
No Matches
CoordinateModel.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//
46// Testing of CoordinateModel
47//
48
49#include "Kokkos_StaticCrsGraph.hpp"
50#include "Kokkos_UnorderedMap.hpp"
54
55#include <set>
56#include <bitset>
57
58#include <Teuchos_Comm.hpp>
59#include <Teuchos_CommHelpers.hpp>
60#include <Teuchos_DefaultComm.hpp>
61#include <Teuchos_ArrayView.hpp>
62#include <Teuchos_OrdinalTraits.hpp>
63
64#include <Tpetra_CrsMatrix.hpp>
65
66using Teuchos::RCP;
67using Teuchos::Comm;
68
69void testCoordinateModel(std::string &fname, int nWeights,
70 const RCP<const Comm<int> > &comm,
71 bool nodeZeroHasAll, bool printInfo)
72{
73 int fail = 0, gfail = 0;
74
75 if (printInfo){
76 std::cout << "Test: " << fname << std::endl;
77 std::cout << "Num Weights: " << nWeights;
78 std::cout << " proc 0 has all: " << nodeZeroHasAll;
79 std::cout << std::endl;
80 }
81
83 // Input data
85
86 typedef Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t> mv_t;
87
88 RCP<UserInputForTests> uinput;
89
90 try{
91 uinput = rcp(new UserInputForTests(testDataFilePath, fname, comm, true));
92 }
93 catch(std::exception &e){
94 fail=1;
95 }
96
97 TEST_FAIL_AND_EXIT(*comm, !fail, "input constructor", 1);
98
99 RCP<mv_t> coords;
100
101 try{
102 coords = uinput->getUICoordinates();
103 }
104 catch(std::exception &e){
105 fail=2;
106 }
107
108 TEST_FAIL_AND_EXIT(*comm, !fail, "getting coordinates", 1);
109
110 int coordDim = coords->getNumVectors();
111
112 TEST_FAIL_AND_EXIT(*comm, coordDim <= 3, "dim 3 at most", 1);
113
114 const zscalar_t *x=NULL, *y=NULL, *z=NULL;
115
116 x = coords->getData(0).getRawPtr();
117 if (coordDim > 1){
118 y = coords->getData(1).getRawPtr();
119 if (coordDim > 2)
120 z = coords->getData(2).getRawPtr();
121 }
122
123 // Are these coordinates correct
124
125 int nLocalIds = coords->getLocalLength();
126 ArrayView<const zgno_t> idList = coords->getMap()->getLocalElementList();
127
128 int nGlobalIds = 0;
129 if (nodeZeroHasAll){
130 if (comm->getRank() > 0){
131 x = y = z = NULL;
132 nLocalIds = 0;
133 }
134 else{
135 nGlobalIds = nLocalIds;
136 }
137 Teuchos::broadcast<int, int>(*comm, 0, &nGlobalIds);
138 }
139 else{
140 nGlobalIds = coords->getGlobalLength();
141 }
142
143 Array<ArrayRCP<const zscalar_t> > coordWeights(nWeights);
144
145 if (nLocalIds > 0){
146 for (int wdim=0; wdim < nWeights; wdim++){
147 zscalar_t *w = new zscalar_t [nLocalIds];
148 for (int i=0; i < nLocalIds; i++){
149 w[i] = ((i%2) + 1) + wdim;
150 }
151 coordWeights[wdim] = Teuchos::arcp(w, 0, nLocalIds);
152 }
153 }
154
155
157 // Create a BasicVectorAdapter adapter object.
159
161 typedef Zoltan2::VectorAdapter<mv_t> base_ia_t;
162
163 RCP<ia_t> ia;
164
165 if (nWeights == 0){ // use the simpler constructor
166 try{
167 ia = rcp(new ia_t(nLocalIds, idList.getRawPtr(), x, y, z));
168 }
169 catch(std::exception &e){
170 fail=3;
171 }
172 }
173 else{
174 std::vector<const zscalar_t *> values, weights;
175 std::vector<int> valueStrides, weightStrides; // default is 1
176 values.push_back(x);
177 if (y) {
178 values.push_back(y);
179 if (z)
180 values.push_back(z);
181 }
182 for (int wdim=0; wdim < nWeights; wdim++){
183 weights.push_back(coordWeights[wdim].getRawPtr());
184 }
185
186 try{
187 ia = rcp(new ia_t(nLocalIds, idList.getRawPtr(),
188 values, valueStrides, weights, weightStrides));
189 }
190 catch(std::exception &e){
191 fail=4;
192 }
193 }
194
195 RCP<const base_ia_t> base_ia = Teuchos::rcp_dynamic_cast<const base_ia_t>(ia);
196
197 TEST_FAIL_AND_EXIT(*comm, !fail, "making input adapter", 1);
198
200 // Create an CoordinateModel with this input
202
204 typedef std::bitset<Zoltan2::NUM_MODEL_FLAGS> modelFlags_t;
206 modelFlags_t modelFlags;
207
208 RCP<const Zoltan2::Environment> env = rcp(new Zoltan2::Environment(comm));
209 RCP<model_t> model;
210
211
212 try{
213 model = rcp(new model_t(base_ia, env, comm, modelFlags));
214 }
215 catch (std::exception &e){
216 fail = 5;
217 }
218
219 TEST_FAIL_AND_EXIT(*comm, !fail, "making model", 1);
220
221 // Test the CoordinateModel interface
222
223 if (model->getCoordinateDim() != coordDim)
224 fail = 6;
225
226 if (!fail && model->getLocalNumCoordinates() != size_t(nLocalIds))
227 fail = 7;
228
229 if (!fail && model->getGlobalNumCoordinates() != size_t(nGlobalIds))
230 fail = 8;
231
232 if (!fail && model->getNumWeightsPerCoordinate() != nWeights)
233 fail = 9;
234
235 gfail = globalFail(*comm, fail);
236
237 if (gfail)
238 printFailureCode(*comm, fail);
239
240 ArrayView<const zgno_t> gids;
241 ArrayView<input_t> xyz;
242 ArrayView<input_t> wgts;
243
244 model->getCoordinates(gids, xyz, wgts);
245
246 if (!fail && gids.size() != nLocalIds)
247 fail = 10;
248
249 for (int i=0; !fail && i < nLocalIds; i++){
250 if (gids[i] != idList[i])
251 fail = 11;
252 }
253
254 if (!fail && wgts.size() != nWeights)
255 fail = 12;
256
257 const zscalar_t *vals[3] = {x, y, z};
258
259 for (int dim=0; !fail && dim < coordDim; dim++){
260 for (int i=0; !fail && i < nLocalIds; i++){
261 if (xyz[dim][i] != vals[dim][i])
262 fail = 13;
263 }
264 }
265
266 for (int wdim=0; !fail && wdim < nWeights; wdim++){
267 for (int i=0; !fail && i < nLocalIds; i++){
268 if (wgts[wdim][i] != coordWeights[wdim][i])
269 fail = 14;
270 }
271 }
272
273 gfail = globalFail(*comm, fail);
274
275 if (gfail)
276 printFailureCode(*comm, fail);
277
278
282 Kokkos::View<const zgno_t *, typename znode_t::device_type> gidsKokkos;
283
284 Kokkos::View<zscalar_t **, Kokkos::LayoutLeft, typename znode_t::device_type> xyzKokkos;
285 Kokkos::View<zscalar_t **, typename znode_t::device_type> wgtsKokkos;
286
287 model->getCoordinatesKokkos(gidsKokkos, xyzKokkos, wgtsKokkos);
288
289 if (!fail && gidsKokkos.extent(0) != static_cast<size_t>(nLocalIds))
290 fail = 10;
291
292 auto gidsKokkos_host = Kokkos::create_mirror_view(gidsKokkos);
293 Kokkos::deep_copy(gidsKokkos_host, gidsKokkos);
294
295 for (int i=0; !fail && i < nLocalIds; i++){
296 if (gidsKokkos_host(i) != idList[i])
297 fail = 11;
298 }
299
300 if (!fail && wgtsKokkos.extent(1) != static_cast<size_t>(nWeights))
301 fail = 12;
302
303 auto xyzKokkos_host = Kokkos::create_mirror_view(xyzKokkos);
304 Kokkos::deep_copy(xyzKokkos_host, xyzKokkos);
305
306 for (int dim=0; !fail && dim < coordDim; dim++){
307 for (int i=0; !fail && i < nLocalIds; i++){
308 if (xyzKokkos_host(i, dim) != vals[dim][i])
309 fail = 13;
310 }
311 }
312
313 auto wgtsKokkos_host = Kokkos::create_mirror_view(wgtsKokkos);
314 Kokkos::deep_copy(wgtsKokkos_host, wgtsKokkos);
315
316 for (int wdim=0; !fail && wdim < nWeights; wdim++){
317 for (int i=0; !fail && i < nLocalIds; i++){
318 if (wgtsKokkos_host(i, wdim) != coordWeights[wdim][i])
319 fail = 14;
320 }
321 }
322
323 gfail = globalFail(*comm, fail);
324
325 if (gfail)
326 printFailureCode(*comm, fail);
327}
328
329int main(int narg, char *arg[])
330{
331 Tpetra::ScopeGuard tscope(&narg, &arg);
332 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
333
334 int rank = comm->getRank();
335 string fname("simple"); // reader will seek coord file
336
337 testCoordinateModel(fname, 0, comm, false, rank==0);
338
339 testCoordinateModel(fname, 1, comm, false, rank==0);
340
341 testCoordinateModel(fname, 2, comm, false, rank==0);
342
343 testCoordinateModel(fname, 0, comm, true, rank==0);
344
345 testCoordinateModel(fname, 1, comm, true, rank==0);
346
347 testCoordinateModel(fname, 2, comm, true, rank==0);
348
349 if (rank==0) std::cout << "PASS" << std::endl;
350
351 return 0;
352}
void testCoordinateModel(std::string &fname, int nWeights, const RCP< const Comm< int > > &comm, bool nodeZeroHasAll, bool printInfo)
int globalFail(const Comm< int > &comm, int fail)
void printFailureCode(const Comm< int > &comm, int fail)
#define TEST_FAIL_AND_EXIT(comm, ok, s, code)
Defines the BasicVectorAdapter class.
Defines the CoordinateModel classes.
common code used by tests
float zscalar_t
std::string testDataFilePath(".")
int main()
BasicVectorAdapter represents a vector (plus optional weights) supplied by the user as pointers to st...
This class provides geometric coordinates with optional weights to the Zoltan2 algorithm.
The user parameters, debug, timing and memory profiling output objects, and error checking methods.
The StridedData class manages lists of weights or coordinates.
VectorAdapter defines the interface for vector input.
static const std::string fail
static ArrayRCP< ArrayRCP< zscalar_t > > weights