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