Zoltan2
BasicCoordinateInput.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 // Test for Zoltan2::BasicVectorAdapter for coordinate-based problems
47 
49 #include <Zoltan2_TestHelpers.hpp>
50 
51 #include <Teuchos_GlobalMPISession.hpp>
52 #include <Teuchos_DefaultComm.hpp>
53 #include <Teuchos_RCP.hpp>
54 #include <Teuchos_CommHelpers.hpp>
55 
56 using Teuchos::RCP;
57 using Teuchos::Comm;
58 using Teuchos::DefaultComm;
59 using Teuchos::Array;
60 
62 
65  int len, int glen, zgno_t *ids,
66  zscalar_t *xyz,
68  int nCoords, int nWeights)
69 {
70  int fail = 0;
71 
72  if (ia->getNumEntriesPerID() != nCoords)
73  fail = 100;
74 
75  if (!fail && ia->getNumWeightsPerID() != nWeights)
76  fail = 101;
77 
78  if (!fail && ia->getLocalNumIDs() != size_t(len))
79  fail = 102;
80 
81  for (int x=0; !fail && x < nCoords; x++){
82  const zgno_t *idList;
83  const zscalar_t *vals;
84  int stride;
85 
86  ia->getIDsView(idList);
87  ia->getEntriesView(vals, stride, x);
88 
89  zscalar_t *coordVal = xyz + x;
90  for (int i=0; !fail && i < len; i++, coordVal += 3){
91 
92  if (idList[i] != ids[i])
93  fail = 105;
94 
95  if (!fail && vals[stride*i] != *coordVal)
96  fail = 106;
97  }
98  }
99 
100  for (int w=0; !fail && w < nWeights; w++){
101  const zscalar_t *wgts;
102  int stride;
103 
104  ia->getWeightsView(wgts, stride, w);
105 
106  zscalar_t *weightVal = weights + len*w;
107  for (int i=0; !fail && i < len; i++, weightVal++){
108  if (wgts[stride*i] != *weightVal)
109  fail = 110;
110  }
111  }
112 
113  return fail;
114 }
115 
116 
117 int main(int argc, char *argv[])
118 {
119  Teuchos::GlobalMPISession session(&argc, &argv);
120  RCP<const Comm<int> > comm = DefaultComm<int>::getComm();
121  int rank = comm->getRank();
122  int fail = 0;
123 
124  // Get some coordinates
125 
126  typedef Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t> mv_t;
127  RCP<UserInputForTests> uinput;
128  std::string fname("simple");
129 
130  try{
131  uinput = rcp(new UserInputForTests(testDataFilePath, fname, comm, true));
132  }
133  catch(std::exception &e){
134  fail=1;
135  }
136 
137  TEST_FAIL_AND_EXIT(*comm, !fail, "input constructor", 1);
138 
139  RCP<mv_t> coords;
140 
141  try{
142  coords = uinput->getUICoordinates();
143  }
144  catch(std::exception &e){
145  fail=1;
146  }
147 
148  TEST_FAIL_AND_EXIT(*comm, !fail, "getting coordinates", 1);
149 
150  int numLocalIds = coords->getLocalLength();
151  int numGlobalIds = coords->getGlobalLength();
152  int coordDim = coords->getNumVectors();
153  ArrayView<const zgno_t> idList = coords->getMap()->getNodeElementList();
154 
155  // Create global Ids, x-, y- and z-coordinates, and also arrays of weights.
156 
157  Array<zgno_t> myIds(numLocalIds);
158  zgno_t base = rank * numLocalIds;
159 
160  int wdim = 2;
161  Array<zscalar_t> weights(numLocalIds*wdim);
162  for (int i = 0; i < numLocalIds*wdim; i++) weights[i] = zscalar_t(i);
163 
164  zscalar_t *x_values= coords->getDataNonConst(0).getRawPtr();
165  zscalar_t *y_values= x_values; // fake 3 dimensions if needed
166  zscalar_t *z_values= x_values;
167 
168  if (coordDim > 1){
169  y_values= coords->getDataNonConst(1).getRawPtr();
170  if (coordDim > 2)
171  z_values= coords->getDataNonConst(2).getRawPtr();
172  }
173 
174  Array<zscalar_t> xyz_values(3*numLocalIds);
175 
176  for (zlno_t i=0; i < numLocalIds; i++) // global Ids
177  myIds[i] = base+i;
178 
179  zscalar_t *x = xyz_values.getRawPtr(); // a stride-3 coordinate array
180  zscalar_t *y = x+1;
181  zscalar_t *z = y+1;
182 
183  for (int i=0, ii=0; i < numLocalIds; i++, ii += 3){
184  x[ii] = x_values[i];
185  y[ii] = y_values[i];
186  z[ii] = z_values[i];
187  }
188 
189  RCP<Zoltan2::BasicVectorAdapter<userTypes_t> > ia;
190 
191  {
193  // 3-dimensional coordinates with stride one and no weights,
194  // using simpler constructor
195 
196  int ncoords = 3;
197  int nweights = 0;
198 
199  try{
201  numLocalIds, myIds.getRawPtr(), x_values, y_values, z_values));
202  }
203  catch (std::exception &e){
204  fail = 1;
205  }
206 
207  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "constructor 0", fail);
208 
209  fail = checkBasicCoordinate(ia.getRawPtr(), numLocalIds, numGlobalIds,
210  myIds.getRawPtr(), xyz_values.getRawPtr(),
211  weights.getRawPtr(), ncoords, nweights);
212 
213  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "check adapter 0", fail);
214  }
215 
216  {
218  // 3-dimensional coordinates with stride one and one weight
219  // using simpler constructor
220 
221  int ncoords = 3;
222  int nweights = 1;
223 
224  try{
226  numLocalIds, myIds.getRawPtr(),
227  x_values, y_values, z_values, 1, 1, 1,
228  true, weights.getRawPtr(), 1));
229  }
230  catch (std::exception &e){
231  fail = 1;
232  }
233 
234  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "constructor 0a", fail);
235 
236  fail = checkBasicCoordinate(ia.getRawPtr(), numLocalIds, numGlobalIds,
237  myIds.getRawPtr(), xyz_values.getRawPtr(),
238  weights.getRawPtr(), ncoords, nweights);
239 
240  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "check adapter 0a", fail);
241  }
242 
243  {
245  // 3-dimensional coordinates with stride one and no weights
246 
247  int ncoords = 3;
248  int nweights = 0;
249 
250  std::vector<const zscalar_t *> values, weightValues;
251  std::vector<int> valueStrides, weightStrides;
252 
253  values.push_back(x_values);
254  values.push_back(y_values);
255  values.push_back(z_values);
256  valueStrides.push_back(1);
257  valueStrides.push_back(1);
258  valueStrides.push_back(1);
259 
260  try{
262  numLocalIds, myIds.getRawPtr(), values, valueStrides,
263  weightValues, weightStrides));
264  }
265  catch (std::exception &e){
266  fail = 1;
267  }
268 
269  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "constructor 1", fail);
270 
271  fail = checkBasicCoordinate(ia.getRawPtr(), numLocalIds, numGlobalIds,
272  myIds.getRawPtr(), xyz_values.getRawPtr(),
273  weights.getRawPtr(), ncoords, nweights);
274 
275  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "check adapter 1", fail);
276 
277  // Try using the default: no strides supplied means strides are one.
278 
279  std::vector<int> emptyStrides;
280 
281  try{
283  numLocalIds, myIds.getRawPtr(), values, emptyStrides,
284  weightValues, emptyStrides));
285  }
286  catch (std::exception &e){
287  fail = 1;
288  }
289 
290  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "constructor 2", fail);
291 
292  fail = checkBasicCoordinate(ia.getRawPtr(), numLocalIds, numGlobalIds,
293  myIds.getRawPtr(), xyz_values.getRawPtr(),
294  weights.getRawPtr(), ncoords, nweights);
295 
296  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "check adapter 2", fail);
297  }
298 
299  {
301  // 2-dimensional coordinates with stride three and two weights
302 
303  int ncoords = 2;
304  int nweights = 2;
305 
306  std::vector<const zscalar_t *> values, weightValues;
307  std::vector<int> valueStrides, weightStrides;
308 
309  values.push_back(xyz_values.getRawPtr());
310  values.push_back(xyz_values.getRawPtr() + 1);
311  valueStrides.push_back(3);
312  valueStrides.push_back(3);
313 
314  weightValues.push_back(weights.getRawPtr());
315  weightValues.push_back(weights.getRawPtr() + numLocalIds);
316  weightStrides.push_back(1);
317  weightStrides.push_back(1);
318 
319  try{
321  numLocalIds, myIds.getRawPtr(), values, valueStrides,
322  weightValues, weightStrides));
323  }
324  catch (std::exception &e){
325  fail = 1;
326  }
327 
328  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "constructor 3", fail);
329 
330  fail = checkBasicCoordinate(ia.getRawPtr(), numLocalIds, numGlobalIds,
331  myIds.getRawPtr(), xyz_values.getRawPtr(),
332  weights.getRawPtr(), ncoords, nweights);
333 
334  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "check adapter 3", fail);
335 
336  // Try using default weight strides
337 
338  std::vector<int> emptyStrides;
339 
340  try{
342  numLocalIds, myIds.getRawPtr(), values, valueStrides,
343  weightValues, emptyStrides));
344  }
345  catch (std::exception &e){
346  fail = 1;
347  }
348 
349  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "constructor 4", fail);
350 
351  fail = checkBasicCoordinate(ia.getRawPtr(), numLocalIds, numGlobalIds,
352  myIds.getRawPtr(), xyz_values.getRawPtr(),
353  weights.getRawPtr(), ncoords, nweights);
354 
355  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "check adapter 4", fail);
356  }
357 
358  {
360  // 1-dimensional coordinates with stride one and two weights
361 
362  int ncoords = 1;
363  int nweights = 2;
364 
365  std::vector<const zscalar_t *> values, weightValues;
366  std::vector<int> valueStrides, weightStrides;
367 
368  values.push_back(x_values);
369  valueStrides.push_back(1);
370 
371  weightValues.push_back(weights.getRawPtr());
372  weightValues.push_back(weights.getRawPtr() + numLocalIds);
373  weightStrides.push_back(1);
374  weightStrides.push_back(1);
375 
376  try{
378  numLocalIds, myIds.getRawPtr(), values, valueStrides,
379  weightValues, weightStrides));
380  }
381  catch (std::exception &e){
382  fail = 1;
383  }
384 
385  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "constructor 4", fail);
386 
387  fail = checkBasicCoordinate(ia.getRawPtr(), numLocalIds, numGlobalIds,
388  myIds.getRawPtr(), xyz_values.getRawPtr(),
389  weights.getRawPtr(), ncoords, nweights);
390 
391  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "check adapter 4", fail);
392  }
393 
394  if (rank == 0)
395  std::cout << "PASS" << std::endl;
396 
397  return fail;
398 }
399 
size_t getLocalNumIDs() const
Returns the number of objects on this process.
int checkBasicCoordinate(Zoltan2::BasicVectorAdapter< userTypes_t > *ia, int len, int glen, zgno_t *ids, zscalar_t *xyz, zscalar_t *weights, int nCoords, int nWeights)
double zscalar_t
A simple class that can be the User template argument for an InputAdapter.
#define TEST_FAIL_AND_EXIT(comm, ok, s, code)
static ArrayRCP< ArrayRCP< zscalar_t > > weights
int zlno_t
common code used by tests
Zoltan2::BasicUserTypes< zscalar_t, zlno_t, zgno_t > userTypes_t
int getNumEntriesPerID() const
Return the number of vectors (typically one).
list idList
Match up parameters to validators.
int main(int argc, char *argv[])
dictionary vals
Definition: xml2dox.py:186
void getWeightsView(const scalar_t *&weights, int &stride, int idx) const
BasicVectorAdapter represents a vector (plus optional weights) supplied by the user as pointers to st...
static const std::string fail
#define TEST_FAIL_AND_RETURN_VALUE(comm, ok, s, rc)
int zgno_t
void getIDsView(const gno_t *&ids) const
int getNumWeightsPerID() const
Returns the number of weights per object. Number of weights per object should be zero or greater...
void getEntriesView(const scalar_t *&entries, int &stride, int idx=0) const
Defines the BasicVectorAdapter class.
std::string testDataFilePath(".")