Zoltan2
XpetraVectorInput.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 
52 #include <string>
53 
55 #include <Zoltan2_InputTraits.hpp>
56 #include <Zoltan2_TestHelpers.hpp>
57 
58 #include <Teuchos_GlobalMPISession.hpp>
59 #include <Teuchos_DefaultComm.hpp>
60 #include <Teuchos_RCP.hpp>
61 #include <Teuchos_Comm.hpp>
62 #include <Teuchos_CommHelpers.hpp>
63 
64 using namespace std;
65 using Teuchos::RCP;
66 using Teuchos::rcp;
67 using Teuchos::rcp_const_cast;
68 using Teuchos::Comm;
69 using Teuchos::DefaultComm;
70 
71 typedef Tpetra::Vector<zscalar_t, zlno_t, zgno_t, znode_t> tvector_t;
72 typedef Xpetra::Vector<zscalar_t, zlno_t, zgno_t, znode_t> xvector_t;
73 
74 void printVector(RCP<const Comm<int> > &comm, zlno_t vlen,
75  const zgno_t *vtxIds, const zscalar_t *vals)
76 {
77  int rank = comm->getRank();
78  int nprocs = comm->getSize();
79  comm->barrier();
80  for (int p=0; p < nprocs; p++){
81  if (p == rank){
82  std::cout << rank << ":" << std::endl;
83  for (zlno_t i=0; i < vlen; i++){
84  std::cout << " " << vtxIds[i] << ": " << vals[i] << std::endl;
85  }
86  std::cout.flush();
87  }
88  comm->barrier();
89  }
90  comm->barrier();
91 }
92 
93 template <typename User>
95  Zoltan2::XpetraMultiVectorAdapter<User> &ia, tvector_t &vector, int wdim,
96  zscalar_t **weights, int *strides)
97 {
98  RCP<const Comm<int> > comm = vector.getMap()->getComm();
99  int fail = 0, gfail=0;
100 
101  if (!fail && ia.getNumEntriesPerID() !=1)
102  fail = 42;
103 
104  if (!fail && ia.getNumWeightsPerID() !=wdim)
105  fail = 41;
106 
107  if (!fail && ia.getLocalNumIDs() != vector.getLocalLength())
108  fail = 4;
109 
110  gfail = globalFail(comm, fail);
111 
112  if (!gfail){
113  const zgno_t *vtxIds=NULL;
114  const zscalar_t *vals=NULL;
115  int stride;
116 
117  size_t nvals = ia.getLocalNumIDs();
118  if (nvals != vector.getLocalLength())
119  fail = 8;
120 
121  ia.getIDsView(vtxIds);
122  ia.getEntriesView(vals, stride);
123 
124  if (!fail && stride != 1)
125  fail = 10;
126 
127  gfail = globalFail(comm, fail);
128 
129  if (gfail == 0){
130  printVector(comm, nvals, vtxIds, vals);
131  }
132  }
133 
134  if (!gfail && wdim){
135  const zscalar_t *wgt =NULL;
136  int stride;
137 
138  for (int w=0; !fail && w < wdim; w++){
139  ia.getWeightsView(wgt, stride, w);
140 
141  if (!fail && stride != strides[w])
142  fail = 101;
143 
144  for (size_t v=0; !fail && v < vector.getLocalLength(); v++){
145  if (wgt[v*stride] != weights[w][v*stride])
146  fail=102;
147  }
148  }
149 
150  gfail = globalFail(comm, fail);
151  }
152 
153  return gfail;
154 }
155 
156 int main(int argc, char *argv[])
157 {
158  Teuchos::GlobalMPISession session(&argc, &argv);
159  RCP<const Comm<int> > comm = DefaultComm<int>::getComm();
160  int rank = comm->getRank();
161  int fail = 0, gfail=0;
162  bool aok = true;
163 
164  // Create object that can give us test Tpetra, Xpetra
165  // and Epetra vectors for testing.
166 
167  RCP<UserInputForTests> uinput;
168 
169  try{
170  uinput =
171  rcp(new UserInputForTests(testDataFilePath,std::string("simple"), comm, true));
172  }
173  catch(std::exception &e){
174  aok = false;
175  std::cout << e.what() << std::endl;
176  }
177  TEST_FAIL_AND_EXIT(*comm, aok, "input ", 1);
178 
179  RCP<tvector_t> tV; // original vector (for checking)
180  RCP<tvector_t> newV; // migrated vector
181 
182  tV = rcp(new tvector_t(uinput->getUITpetraCrsGraph()->getRowMap()));
183  tV->randomize();
184  size_t vlen = tV->getLocalLength();
185 
186  // To test migration in the input adapter we need a Solution object.
187 
188  RCP<const Zoltan2::Environment> env = rcp(new Zoltan2::Environment);
189 
190  int nWeights = 1;
191 
193  typedef adapter_t::part_t part_t;
194 
195  part_t *p = new part_t [vlen];
196  memset(p, 0, sizeof(part_t) * vlen);
197  ArrayRCP<part_t> solnParts(p, 0, vlen, true);
198 
199  std::vector<const zscalar_t *> emptyWeights;
200  std::vector<int> emptyStrides;
201 
202  Zoltan2::PartitioningSolution<adapter_t> solution(env, comm, nWeights);
203  solution.setParts(solnParts);
204 
206  // User object is Tpetra::Vector, no weights
207  if (!gfail){
208  if (rank==0)
209  std::cout << "Constructed with Tpetra::Vector" << std::endl;
210 
211  RCP<const tvector_t> ctV = rcp_const_cast<const tvector_t>(tV);
212  RCP<adapter_t> tVInput;
213 
214  try {
215  tVInput =
217  emptyWeights, emptyStrides));
218  }
219  catch (std::exception &e){
220  aok = false;
221  std::cout << e.what() << std::endl;
222  }
223  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter ", 1);
224 
225  fail = verifyInputAdapter<tvector_t>(*tVInput, *tV, 0, NULL, NULL);
226 
227  gfail = globalFail(comm, fail);
228 
229  if (!gfail){
230  tvector_t *vMigrate = NULL;
231  try{
232  tVInput->applyPartitioningSolution(*tV, vMigrate, solution);
233  newV = rcp(vMigrate);
234  }
235  catch (std::exception &e){
236  fail = 11;
237  }
238 
239  gfail = globalFail(comm, fail);
240 
241  if (!gfail){
242  RCP<const tvector_t> cnewV = rcp_const_cast<const tvector_t>(newV);
243  RCP<Zoltan2::XpetraMultiVectorAdapter<tvector_t> > newInput;
244  try{
245  newInput = rcp(new Zoltan2::XpetraMultiVectorAdapter<tvector_t>(cnewV,
246  emptyWeights, emptyStrides));
247  }
248  catch (std::exception &e){
249  aok = false;
250  std::cout << e.what() << std::endl;
251  }
252  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 2 ", 1);
253 
254  if (rank==0){
255  std::cout << "Constructed with ";
256  std::cout << "Tpetra::Vector migrated to proc 0" << std::endl;
257  }
258  fail = verifyInputAdapter<tvector_t>(*newInput, *newV, 0, NULL, NULL);
259  if (fail) fail += 100;
260  gfail = globalFail(comm, fail);
261  }
262  }
263  if (gfail){
264  printFailureCode(comm, fail);
265  }
266  }
267 
269  // User object is Xpetra::Vector
270  if (!gfail){
271  if (rank==0)
272  std::cout << "Constructed with Xpetra::Vector" << std::endl;
273 
274  RCP<tvector_t> xtV =
275  rcp(new tvector_t(uinput->getUITpetraCrsGraph()->getRowMap()));
276  xtV->randomize();
277  RCP<xvector_t> xV = Zoltan2::XpetraTraits<tvector_t>::convertToXpetra(xtV);
278  RCP<const xvector_t> cxV = rcp_const_cast<const xvector_t>(xV);
279  RCP<Zoltan2::XpetraMultiVectorAdapter<xvector_t> > xVInput;
280 
281  try {
282  xVInput =
284  emptyWeights, emptyStrides));
285  }
286  catch (std::exception &e){
287  aok = false;
288  std::cout << e.what() << std::endl;
289  }
290  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 3 ", 1);
291 
292  fail = verifyInputAdapter<xvector_t>(*xVInput, *tV, 0, NULL, NULL);
293 
294  gfail = globalFail(comm, fail);
295 
296  if (!gfail){
297  xvector_t *vMigrate =NULL;
298  try{
299  xVInput->applyPartitioningSolution(*xV, vMigrate, solution);
300  }
301  catch (std::exception &e){
302  fail = 11;
303  }
304 
305  gfail = globalFail(comm, fail);
306 
307  if (!gfail){
308  RCP<const xvector_t> cnewV(vMigrate);
309  RCP<Zoltan2::XpetraMultiVectorAdapter<xvector_t> > newInput;
310  try{
311  newInput =
313  emptyWeights, emptyStrides));
314  }
315  catch (std::exception &e){
316  aok = false;
317  std::cout << e.what() << std::endl;
318  }
319  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 4 ", 1);
320 
321  if (rank==0){
322  std::cout << "Constructed with ";
323  std::cout << "Xpetra::Vector migrated to proc 0" << std::endl;
324  }
325  fail = verifyInputAdapter<xvector_t>(*newInput, *newV, 0, NULL, NULL);
326  if (fail) fail += 100;
327  gfail = globalFail(comm, fail);
328  }
329  }
330  if (gfail){
331  printFailureCode(comm, fail);
332  }
333  }
334 
335 #ifdef HAVE_EPETRA_DATA_TYPES
336  // User object is Epetra_Vector
338  typedef Epetra_Vector evector_t;
339  if (!gfail){
340  if (rank==0)
341  std::cout << "Constructed with Epetra_Vector" << std::endl;
342 
343  RCP<evector_t> eV =
344  rcp(new Epetra_Vector(uinput->getUIEpetraCrsGraph()->RowMap()));
345  eV->Random();
346  RCP<const evector_t> ceV = rcp_const_cast<const evector_t>(eV);
347  RCP<Zoltan2::XpetraMultiVectorAdapter<evector_t> > eVInput;
348 
349  bool goodAdapter = true;
350  try {
351  eVInput =
353  emptyWeights, emptyStrides));
354  }
355  catch (std::exception &e){
356  if (std::is_same<znode_t, Xpetra::EpetraNode>::value) {
357  aok = false;
358  goodAdapter = false;
359  std::cout << e.what() << std::endl;
360  }
361  else {
362  // We expect an error from Xpetra when znode_t != Xpetra::EpetraNode
363  // Ignore it, but skip tests using this matrix adapter.
364  std::cout << "Node type is not supported by Xpetra's Epetra interface;"
365  << " Skipping this test." << std::endl;
366  std::cout << "FYI: Here's the exception message: " << std::endl
367  << e.what() << std::endl;
368  goodAdapter = false;
369  }
370  }
371  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 5 ", 1);
372 
373  if (goodAdapter) {
374  fail = verifyInputAdapter<evector_t>(*eVInput, *tV, 0, NULL, NULL);
375 
376  gfail = globalFail(comm, fail);
377 
378  if (!gfail){
379  evector_t *vMigrate =NULL;
380  try{
381  eVInput->applyPartitioningSolution(*eV, vMigrate, solution);
382  }
383  catch (std::exception &e){
384  fail = 11;
385  }
386 
387  gfail = globalFail(comm, fail);
388 
389  if (!gfail){
390  RCP<const evector_t> cnewV(vMigrate, true);
391  RCP<Zoltan2::XpetraMultiVectorAdapter<evector_t> > newInput;
392  try{
393  newInput =
395  emptyWeights, emptyStrides));
396  }
397  catch (std::exception &e){
398  aok = false;
399  std::cout << e.what() << std::endl;
400  }
401  TEST_FAIL_AND_EXIT(*comm, aok, "XpetraMultiVectorAdapter 6 ", 1);
402 
403  if (rank==0){
404  std::cout << "Constructed with ";
405  std::cout << "Epetra_Vector migrated to proc 0" << std::endl;
406  }
407  fail = verifyInputAdapter<evector_t>(*newInput, *newV, 0, NULL, NULL);
408  if (fail) fail += 100;
409  gfail = globalFail(comm, fail);
410  }
411  }
412  if (gfail){
413  printFailureCode(comm, fail);
414  }
415  }
416  }
417 #endif
418 
420  // DONE
421 
422  if (rank==0)
423  std::cout << "PASS" << std::endl;
424 }
int globalFail(const RCP< const Comm< int > > &comm, int fail)
double zscalar_t
int main(int argc, char *argv[])
#define TEST_FAIL_AND_EXIT(comm, ok, s, code)
Tpetra::Vector< zscalar_t, zlno_t, zgno_t, znode_t > tvector_t
static ArrayRCP< ArrayRCP< zscalar_t > > weights
Tpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > tvector_t
int zlno_t
common code used by tests
static RCP< User > convertToXpetra(const RCP< User > &a)
Convert the object to its Xpetra wrapped version.
Defines the XpetraMultiVectorAdapter.
size_t getLocalNumIDs() const
Returns the number of objects on this process.
int verifyInputAdapter(Zoltan2::XpetraMultiVectorAdapter< User > &ia, tvector_t &vector, int wdim, zscalar_t **weights, int *strides)
dictionary vals
Definition: xml2dox.py:186
Xpetra::Vector< zscalar_t, zlno_t, zgno_t, znode_t > xvector_t
A PartitioningSolution is a solution to a partitioning problem.
int getNumEntriesPerID() const
Return the number of vectors (typically one).
void getEntriesView(const scalar_t *&elements, int &stride, int idx=0) const
void printVector(RCP< const Comm< int > > &comm, zlno_t vlen, const zgno_t *vtxIds, const zscalar_t *vals)
Traits for application input objects.
The user parameters, debug, timing and memory profiling output objects, and error checking methods...
An adapter for Xpetra::MultiVector.
static const std::string fail
void setParts(ArrayRCP< part_t > &partList)
The algorithm uses setParts to set the solution.
int zgno_t
void printFailureCode(const RCP< const Comm< int > > &comm, int fail)
void getIDsView(const gno_t *&ids) const
Provide a pointer to this process&#39; identifiers.
Xpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > xvector_t
int getNumWeightsPerID() const
Returns the number of weights per object. Number of weights per object should be zero or greater...
void getWeightsView(const scalar_t *&weights, int &stride, int idx) const
std::string testDataFilePath(".")