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