Zoltan2
BasicVectorInput.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
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 
61 
63  Zoltan2::BasicVectorAdapter<userTypes_t> *ia, int len, int glen,
64  zgno_t *ids, int mvdim, const zscalar_t **values, int *valueStrides,
65  int wdim, const zscalar_t **weights, int *weightStrides)
66 {
67  int fail = 0;
68  bool strideOne = false;
69 
70  if (valueStrides == NULL) strideOne = true;
71 
72  if (ia->getNumEntriesPerID() != mvdim)
73  fail = 100;
74 
75  if (!fail && ia->getNumWeightsPerID() != wdim)
76  fail = 101;
77 
78  if (!fail && ia->getLocalNumIDs() != size_t(len))
79  fail = 102;
80 
81  const zgno_t *idList;
82  ia->getIDsView(idList);
83  for (int i=0; !fail && i < len; i++)
84  if (!fail && idList[i] != ids[i])
85  fail = 107;
86 
87  for (int v=0; !fail && v < mvdim; v++){
88  const zscalar_t *vals;
89  int correctStride = (strideOne ? 1 : valueStrides[v]);
90  int stride;
91 
92  ia->getEntriesView(vals, stride, v);
93 
94  if (!fail && stride != correctStride)
95  fail = 105;
96 
97  for (int i=0; !fail && i < len; i++){
98 // TODO fix values check
99 // if (vals[stride*i] != values[v][correctStride*i])
100 // fail = 106;
101  }
102  }
103 
104  for (int w=0; !fail && w < wdim; w++){
105  const zscalar_t *wgts;
106  int stride;
107 
108  ia->getWeightsView(wgts, stride, w);
109 
110  if (!fail && stride != weightStrides[w])
111  fail = 109;
112 
113  for (int i=0; !fail && i < len; i++){
114  if (wgts[stride*i] != weights[w][weightStrides[w]*i])
115  fail = 110;
116  }
117  }
118 
119  return fail;
120 }
121 
122 
123 int main(int argc, char *argv[])
124 {
125  Teuchos::GlobalMPISession session(&argc, &argv);
126  RCP<const Comm<int> > comm = DefaultComm<int>::getComm();
127  int rank = comm->getRank();
128  int nprocs = comm->getSize();
129  int fail = 0;
130 
131  // Create a single vector and a strided multi-vector with
132  // strided multi-weights.
133 
134  zlno_t numLocalIds = 10;
135  zgno_t *myIds = new zgno_t[numLocalIds];
136  zgno_t base = rank * numLocalIds;
137 
138  int wdim = 2;
139  zscalar_t *weights = new zscalar_t [numLocalIds*wdim];
140  int *weightStrides = new int [wdim];
141  const zscalar_t **weightPtrs = new const zscalar_t * [wdim];
142 
143  int mvdim = 3;
144  zscalar_t *v_values= new zscalar_t [numLocalIds];
145  zscalar_t *mv_values= new zscalar_t [mvdim * numLocalIds];
146  int *valueStrides = new int [mvdim];
147  const zscalar_t **valuePtrs = new const zscalar_t * [mvdim];
148 
149  for (zlno_t i=0; i < numLocalIds; i++){
150  myIds[i] = base+i;
151 
152  for (int w=0; w < wdim; w++)
153  weights[w*numLocalIds + i] = w + 1 + nprocs - rank;
154 
155  v_values[i] = numLocalIds-i;
156 
157  for (int v=0; v < mvdim; v++)
158  mv_values[i*mvdim + v] = (v+1) * (nprocs-rank) / (i+1);
159  }
160 
161  for (int w=0; w < wdim; w++){
162  weightStrides[w] = 1;
163  weightPtrs[w] = weights + numLocalIds*w;
164  }
165 
166  for (int v=0; v < mvdim; v++){
167  valueStrides[v] = mvdim;
168  valuePtrs[v] = mv_values + v;
169  }
170 
172 
173  {
174  // A Zoltan2::BasicVectorAdapter object with one vector and no weights
175 
176  std::vector<const zscalar_t *> weightValues;
177  std::vector<int> strides;
178 
179  try{
180  ia = new Zoltan2::BasicVectorAdapter<userTypes_t>(numLocalIds, myIds,
181  v_values, 1);
182  }
183  catch (std::exception &e){
184  fail = 1;
185  }
186 
187  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "constructor 1", fail);
188 
189  fail = checkBasicVector(ia, numLocalIds, numLocalIds*nprocs,
190  myIds, 1, valuePtrs, NULL, 0, NULL, NULL);
191 
192  delete ia;
193 
194  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "check vector 1", fail);
195  }
196 
197  {
198  // A Zoltan2::BasicVectorAdapter object with one vector and weights
199 
200  std::vector<const zscalar_t *> weightValues;
201  std::vector<int> strides;
202 
203  weightValues.push_back(weightPtrs[0]);
204  weightValues.push_back(weightPtrs[1]);
205  strides.push_back(1);
206  strides.push_back(1);
207 
208  try{
209  ia = new Zoltan2::BasicVectorAdapter<userTypes_t>(numLocalIds, myIds,
210  v_values, 1, true, weightPtrs[0], 1);
211  }
212  catch (std::exception &e){
213  fail = 1;
214  }
215 
216  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "constructor 2", fail);
217 
218  fail = checkBasicVector(ia, numLocalIds, numLocalIds*nprocs,
219  myIds, 1, valuePtrs, NULL, 1, weightPtrs, weightStrides);
220 
221  delete ia;
222 
223  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "check vector 2", fail);
224  }
225 
226  {
227  // A Zoltan2::BasicVectorAdapter object with a multivector and no weights
228 
229  std::vector<const zscalar_t *> weightValues, values;
230  std::vector<int> wstrides, vstrides;
231 
232  for (int dim=0; dim < mvdim; dim++){
233  values.push_back(valuePtrs[dim]);
234  vstrides.push_back(mvdim);
235  }
236 
237 
238  try{
240  numLocalIds, myIds, values, vstrides, weightValues, wstrides);
241  }
242  catch (std::exception &e){
243  fail = 1;
244  }
245 
246  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "constructor 3", fail);
247 
248  fail = checkBasicVector(ia, numLocalIds, numLocalIds*nprocs,
249  myIds, mvdim, valuePtrs, valueStrides, 0, NULL, NULL);
250 
251  delete ia;
252 
253  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "check vector 3", fail);
254  }
255 
256  {
257  // A Zoltan2::BasicVectorAdapter object with a multivector with weights
258 
259  std::vector<const zscalar_t *> weightValues, values;
260  std::vector<int> wstrides, vstrides;
261 
262  for (int dim=0; dim < wdim; dim++){
263  weightValues.push_back(weightPtrs[dim]);
264  wstrides.push_back(1);
265  }
266 
267  for (int dim=0; dim < mvdim; dim++){
268  values.push_back(valuePtrs[dim]);
269  vstrides.push_back(mvdim);
270  }
271 
272  try{
274  numLocalIds, myIds, values, vstrides, weightValues, wstrides);
275 
276  }
277  catch (std::exception &e){
278  fail = 1;
279  }
280 
281  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "constructor 4", fail);
282 
283  fail = checkBasicVector(ia, numLocalIds, numLocalIds*nprocs,
284  myIds, mvdim, valuePtrs, valueStrides,
285  wdim, weightPtrs, weightStrides);
286 
287  delete ia;
288 
289  TEST_FAIL_AND_RETURN_VALUE(*comm, fail==0, "check vector 4", fail);
290  }
291 
292  if (rank == 0)
293  std::cout << "PASS" << std::endl;
294 
295  delete [] myIds;
296  delete [] weights;
297  delete [] weightStrides;
298  delete [] weightPtrs;
299  delete [] v_values;
300  delete [] mv_values;
301  delete [] valueStrides;
302  delete [] valuePtrs;
303 
304  return fail;
305 }
306 
size_t getLocalNumIDs() const
Returns the number of objects on this process.
double zscalar_t
A simple class that can be the User template argument for an InputAdapter.
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.
dictionary vals
Definition: xml2dox.py:186
void getWeightsView(const scalar_t *&weights, int &stride, int idx) const
int main(int argc, char *argv[])
int checkBasicVector(Zoltan2::BasicVectorAdapter< userTypes_t > *ia, int len, int glen, zgno_t *ids, int mvdim, const zscalar_t **values, int *valueStrides, int wdim, const zscalar_t **weights, int *weightStrides)
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.