Zoltan2
XpetraTraits.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 // Basic test of the XpetraTraits definitions.
47 //
48 // TODO - a real test would figure out if the migrated objects are
49 // the same as the original, here we just look at them on stdout.
50 // TODO look at number of diagonals and max number of entries in
51 // Tpetra and Xpetra migrated graphs. They're garbage.
52 
53 #include <Zoltan2_XpetraTraits.hpp>
54 #include <Zoltan2_TestHelpers.hpp>
55 
56 #include <string>
57 #include <Teuchos_GlobalMPISession.hpp>
58 #include <Teuchos_DefaultComm.hpp>
59 #include <Teuchos_RCP.hpp>
60 #include <Teuchos_Array.hpp>
61 #include <Teuchos_ArrayRCP.hpp>
62 #include <Teuchos_Comm.hpp>
63 #include <Teuchos_VerboseObject.hpp>
64 #include <Tpetra_CrsMatrix.hpp>
65 #include <Tpetra_Vector.hpp>
66 
67 #ifdef HAVE_ZOLTAN_EPETRA
68 #include <Xpetra_EpetraUtils.hpp>
69 #endif
70 
71 using namespace std;
72 using std::string;
73 using Teuchos::RCP;
74 using Teuchos::ArrayRCP;
75 using Teuchos::ArrayView;
76 using Teuchos::Array;
77 using Teuchos::rcp;
78 using Teuchos::Comm;
79 
80 ArrayRCP<zgno_t> roundRobinMap(
81  const RCP<const Xpetra::Map<zlno_t, zgno_t, znode_t> > &m)
82 {
83  const RCP<const Comm<int> > &comm = m->getComm();
84  int proc = comm->getRank();
85  int nprocs = comm->getSize();
86  zgno_t base = m->getMinAllGlobalIndex();
87  zgno_t max = m->getMaxAllGlobalIndex();
88  size_t globalrows = m->getGlobalNumElements();
89  if (globalrows != size_t(max - base + 1)){
90  TEST_FAIL_AND_EXIT(*comm, 0, "Map is invalid for test - fix test", 1);
91  }
92  RCP<Array<zgno_t> > mygids = rcp(new Array<zgno_t>);
93  zgno_t firstzgno_t = proc;
94  if (firstzgno_t < base){
95  zgno_t n = base % proc;
96  if (n>0)
97  firstzgno_t = base - n + proc;
98  else
99  firstzgno_t = base;
100  }
101  for (zgno_t gid=firstzgno_t; gid <= max; gid+=nprocs){
102  (*mygids).append(gid);
103  }
104 
105  ArrayRCP<zgno_t> newIdArcp = Teuchos::arcp(mygids);
106 
107  return newIdArcp;
108 }
109 
110 int main(int argc, char *argv[])
111 {
112  Teuchos::GlobalMPISession session(&argc, &argv);
113  RCP<const Comm<int> > comm = Teuchos::DefaultComm<int>::getComm();
114  int rank = comm->getRank();
115  bool aok = true;
116 
117  Teuchos::RCP<Teuchos::FancyOStream> outStream =
118  Teuchos::VerboseObjectBase::getDefaultOStream();
119  Teuchos::EVerbosityLevel v=Teuchos::VERB_EXTREME;
120 
121  typedef Tpetra::CrsMatrix<zscalar_t,zlno_t,zgno_t,znode_t> tmatrix_t;
122  typedef Tpetra::CrsGraph<zlno_t,zgno_t,znode_t> tgraph_t;
123  typedef Tpetra::Vector<zscalar_t,zlno_t,zgno_t,znode_t> tvector_t;
124  typedef Tpetra::MultiVector<zscalar_t,zlno_t,zgno_t,znode_t> tmvector_t;
125  typedef Xpetra::CrsMatrix<zscalar_t,zlno_t,zgno_t,znode_t> xmatrix_t;
126  typedef Xpetra::CrsGraph<zlno_t,zgno_t,znode_t> xgraph_t;
127  typedef Xpetra::Vector<zscalar_t,zlno_t,zgno_t,znode_t> xvector_t;
128  typedef Xpetra::MultiVector<zscalar_t,zlno_t,zgno_t,znode_t> xmvector_t;
129  typedef Xpetra::TpetraMap<zlno_t,zgno_t,znode_t> xtmap_t;
130 
131  // Create object that can give us test Tpetra and Xpetra input.
132 
133  RCP<UserInputForTests> uinput;
134 
135  try{
136  uinput =
137  rcp(new UserInputForTests(testDataFilePath,std::string("simple"), comm, true));
138  }
139  catch(std::exception &e){
140  aok = false;
141  std::cout << e.what() << std::endl;
142  }
143  TEST_FAIL_AND_EXIT(*comm, aok, "input ", 1);
144 
146  // Tpetra::CrsMatrix
147  // Tpetra::CrsGraph
148  // Tpetra::Vector
149  // Tpetra::MultiVector
151 
152  // XpetraTraits<Tpetra::CrsMatrix<zscalar_t, zlno_t, zgno_t, znode_t> >
153  {
154  RCP<tmatrix_t> M;
155 
156  try{
157  M = uinput->getUITpetraCrsMatrix();
158  }
159  catch(std::exception &e){
160  aok = false;
161  std::cout << e.what() << std::endl;
162  }
163  TEST_FAIL_AND_EXIT(*comm, aok, "getTpetraCrsMatrix ", 1);
164 
165  if (rank== 0)
166  std::cout << "Original Tpetra matrix " << M->getGlobalNumRows()
167  << " x " << M->getGlobalNumCols() << std::endl;
168 
169  M->describe(*outStream,v);
170 
171  RCP<const xtmap_t> xmap(new xtmap_t(M->getRowMap()));
172 
173  ArrayRCP<zgno_t> newRowIds = roundRobinMap(xmap);
174 
175  zgno_t localNumRows = newRowIds.size();
176 
177  RCP<const tmatrix_t> newM;
178  try{
180  localNumRows, newRowIds.getRawPtr());
181  }
182  catch(std::exception &e){
183  aok = false;
184  std::cout << e.what() << std::endl;
185  }
186  TEST_FAIL_AND_EXIT(*comm, aok,
187  " Zoltan2::XpetraTraits<tmatrix_t>::doMigration ", 1);
188 
189  if (rank== 0)
190  std::cout << "Migrated Tpetra matrix" << std::endl;
191 
192  newM->describe(*outStream,v);
193  }
194 
195  // XpetraTraits<Tpetra::CrsGraph<zscalar_t, zlno_t, zgno_t, znode_t> >
196  {
197  RCP<tgraph_t> G;
198 
199  try{
200  G = uinput->getUITpetraCrsGraph();
201  }
202  catch(std::exception &e){
203  aok = false;
204  std::cout << e.what() << std::endl;
205  }
206  TEST_FAIL_AND_EXIT(*comm, aok, "getTpetraCrsGraph ", 1);
207 
208  if (rank== 0)
209  std::cout << "Original Tpetra graph" << std::endl;
210 
211  G->describe(*outStream,v);
212 
213  RCP<const xtmap_t> xmap(new xtmap_t(G->getRowMap()));
214  ArrayRCP<zgno_t> newRowIds = roundRobinMap(xmap);
215 
216  zgno_t localNumRows = newRowIds.size();
217 
218  RCP<const tgraph_t> newG;
219  try{
221  localNumRows, newRowIds.getRawPtr());
222  }
223  catch(std::exception &e){
224  aok = false;
225  std::cout << e.what() << std::endl;
226  }
227  TEST_FAIL_AND_EXIT(*comm, aok,
228  " Zoltan2::XpetraTraits<tgraph_t>::doMigration ", 1);
229 
230  if (rank== 0)
231  std::cout << "Migrated Tpetra graph" << std::endl;
232 
233  newG->describe(*outStream,v);
234  }
235 
236  // XpetraTraits<Tpetra::Vector<zscalar_t, zlno_t, zgno_t, znode_t>>
237  {
238  RCP<tvector_t> V;
239 
240  try{
241  V = rcp(new tvector_t(uinput->getUITpetraCrsGraph()->getRowMap(), 1));
242  V->randomize();
243  }
244  catch(std::exception &e){
245  aok = false;
246  std::cout << e.what() << std::endl;
247  }
248  TEST_FAIL_AND_EXIT(*comm, aok, "getTpetraVector", 1);
249 
250  if (rank== 0)
251  std::cout << "Original Tpetra vector" << std::endl;
252 
253  V->describe(*outStream,v);
254 
255  RCP<const xtmap_t> xmap(new xtmap_t(V->getMap()));
256  ArrayRCP<zgno_t> newRowIds = roundRobinMap(xmap);
257 
258  zgno_t localNumRows = newRowIds.size();
259 
260  RCP<const tvector_t> newV;
261  try{
263  localNumRows, newRowIds.getRawPtr());
264  }
265  catch(std::exception &e){
266  aok = false;
267  std::cout << e.what() << std::endl;
268  }
269  TEST_FAIL_AND_EXIT(*comm, aok,
270  " Zoltan2::XpetraTraits<tvector_t>::doMigration ", 1);
271 
272  if (rank== 0)
273  std::cout << "Migrated Tpetra vector" << std::endl;
274 
275  newV->describe(*outStream,v);
276  }
277 
278  // XpetraTraits<Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t>>
279  {
280  RCP<tmvector_t> MV;
281 
282  try{
283  MV = rcp(new tmvector_t(uinput->getUITpetraCrsGraph()->getRowMap(), 3));
284  MV->randomize();
285  }
286  catch(std::exception &e){
287  aok = false;
288  std::cout << e.what() << std::endl;
289  }
290  TEST_FAIL_AND_EXIT(*comm, aok, "getTpetraMultiVector", 1);
291 
292  if (rank== 0)
293  std::cout << "Original Tpetra multivector" << std::endl;
294 
295  MV->describe(*outStream,v);
296 
297  RCP<const xtmap_t> xmap(new xtmap_t(MV->getMap()));
298  ArrayRCP<zgno_t> newRowIds = roundRobinMap(xmap);
299 
300  zgno_t localNumRows = newRowIds.size();
301 
302  RCP<const tmvector_t> newMV;
303  try{
305  localNumRows, newRowIds.getRawPtr());
306  }
307  catch(std::exception &e){
308  aok = false;
309  std::cout << e.what() << std::endl;
310  }
311  TEST_FAIL_AND_EXIT(*comm, aok,
312  " Zoltan2::XpetraTraits<tmvector_t>::doMigration ", 1);
313 
314  if (rank== 0)
315  std::cout << "Migrated Tpetra multivector" << std::endl;
316 
317  newMV->describe(*outStream,v);
318  }
319 
321  // Xpetra::CrsMatrix
322  // Xpetra::CrsGraph
323  // Xpetra::Vector
324  // Xpetra::MultiVector
326 
327  // XpetraTraits<Xpetra::CrsMatrix<zscalar_t, zlno_t, zgno_t, znode_t> >
328  {
329  RCP<xmatrix_t> M;
330 
331  try{
332  M = uinput->getUIXpetraCrsMatrix();
333  }
334  catch(std::exception &e){
335  aok = false;
336  std::cout << e.what() << std::endl;
337  }
338  TEST_FAIL_AND_EXIT(*comm, aok, "getXpetraCrsMatrix ", 1);
339 
340  if (rank== 0)
341  std::cout << "Original Xpetra matrix" << std::endl;
342 
343  M->describe(*outStream,v);
344 
345  ArrayRCP<zgno_t> newRowIds = roundRobinMap(M->getRowMap());
346 
347  zgno_t localNumRows = newRowIds.size();
348 
349  RCP<const xmatrix_t> newM;
350  try{
352  localNumRows, newRowIds.getRawPtr());
353  }
354  catch(std::exception &e){
355  aok = false;
356  std::cout << e.what() << std::endl;
357  }
358  TEST_FAIL_AND_EXIT(*comm, aok,
359  " Zoltan2::XpetraTraits<xmatrix_t>::doMigration ", 1);
360 
361  if (rank== 0)
362  std::cout << "Migrated Xpetra matrix" << std::endl;
363 
364  newM->describe(*outStream,v);
365  }
366 
367  // XpetraTraits<Xpetra::CrsGraph<zscalar_t, zlno_t, zgno_t, znode_t> >
368  {
369  RCP<xgraph_t> G;
370 
371  try{
372  G = uinput->getUIXpetraCrsGraph();
373  }
374  catch(std::exception &e){
375  aok = false;
376  std::cout << e.what() << std::endl;
377  }
378  TEST_FAIL_AND_EXIT(*comm, aok, "getXpetraCrsGraph ", 1);
379 
380  if (rank== 0)
381  std::cout << "Original Xpetra graph" << std::endl;
382 
383  G->describe(*outStream,v);
384 
385  ArrayRCP<zgno_t> newRowIds = roundRobinMap(G->getRowMap());
386 
387  zgno_t localNumRows = newRowIds.size();
388 
389  RCP<const xgraph_t> newG;
390  try{
392  localNumRows, newRowIds.getRawPtr());
393  }
394  catch(std::exception &e){
395  aok = false;
396  std::cout << e.what() << std::endl;
397  }
398  TEST_FAIL_AND_EXIT(*comm, aok,
399  " Zoltan2::XpetraTraits<xgraph_t>::doMigration ", 1);
400 
401  if (rank== 0)
402  std::cout << "Migrated Xpetra graph" << std::endl;
403 
404  newG->describe(*outStream,v);
405  }
406 
407  // XpetraTraits<Xpetra::Vector<zscalar_t, zlno_t, zgno_t, znode_t>>
408  {
409  RCP<xvector_t> V;
410 
411  try{
412  RCP<tvector_t> tV =
413  rcp(new tvector_t(uinput->getUITpetraCrsGraph()->getRowMap(), 1));
414  tV->randomize();
416  }
417  catch(std::exception &e){
418  aok = false;
419  std::cout << e.what() << std::endl;
420  }
421  TEST_FAIL_AND_EXIT(*comm, aok, "getXpetraVector", 1);
422 
423  if (rank== 0)
424  std::cout << "Original Xpetra vector" << std::endl;
425 
426  V->describe(*outStream,v);
427 
428  ArrayRCP<zgno_t> newRowIds = roundRobinMap(V->getMap());
429 
430  zgno_t localNumRows = newRowIds.size();
431 
432  RCP<const xvector_t> newV;
433  try{
435  localNumRows, newRowIds.getRawPtr());
436  }
437  catch(std::exception &e){
438  aok = false;
439  std::cout << e.what() << std::endl;
440  }
441  TEST_FAIL_AND_EXIT(*comm, aok,
442  " Zoltan2::XpetraTraits<xvector_t>::doMigration ", 1);
443 
444  if (rank== 0)
445  std::cout << "Migrated Xpetra vector" << std::endl;
446 
447  newV->describe(*outStream,v);
448  }
449 
450  // XpetraTraits<Xpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t>>
451  {
452  RCP<xmvector_t> MV;
453 
454  try{
455  RCP<tmvector_t> tMV =
456  rcp(new tmvector_t(uinput->getUITpetraCrsGraph()->getRowMap(), 3));
457  tMV->randomize();
459  }
460  catch(std::exception &e){
461  aok = false;
462  std::cout << e.what() << std::endl;
463  }
464  TEST_FAIL_AND_EXIT(*comm, aok, "getXpetraMultiVector", 1);
465 
466  if (rank== 0)
467  std::cout << "Original Xpetra multivector" << std::endl;
468 
469  MV->describe(*outStream,v);
470 
471  ArrayRCP<zgno_t> newRowIds = roundRobinMap(MV->getMap());
472 
473  zgno_t localNumRows = newRowIds.size();
474 
475  RCP<const xmvector_t> newMV;
476  try{
478  localNumRows, newRowIds.getRawPtr());
479  }
480  catch(std::exception &e){
481  aok = false;
482  std::cout << e.what() << std::endl;
483  }
484  TEST_FAIL_AND_EXIT(*comm, aok,
485  " Zoltan2::XpetraTraits<xmvector_t>::doMigration ", 1);
486 
487  if (rank== 0)
488  std::cout << "Migrated Xpetra multivector" << std::endl;
489 
490  newMV->describe(*outStream,v);
491  }
492 
493 #ifdef HAVE_EPETRA_DATA_TYPES
494  // Epetra_CrsMatrix
496  // Epetra_CrsGraph
497  // Epetra_Vector
498  // Epetra_MultiVector
500 
501  typedef Epetra_CrsMatrix ematrix_t;
502  typedef Epetra_CrsGraph egraph_t;
503  typedef Epetra_Vector evector_t;
504  typedef Epetra_MultiVector emvector_t;
505  typedef Xpetra::EpetraMapT<zgno_t, znode_t> xemap_t;
506  typedef Epetra_BlockMap emap_t;
507 
508  // Create object that can give us test Epetra input.
509 
510  RCP<UserInputForTests> euinput;
511 
512  try{
513  euinput =
514  rcp(new UserInputForTests(testDataFilePath,std::string("simple"), comm, true));
515  }
516  catch(std::exception &e){
517  aok = false;
518  std::cout << e.what() << std::endl;
519  }
520  TEST_FAIL_AND_EXIT(*comm, aok, "epetra input ", 1);
521 
522  // XpetraTraits<Epetra_CrsMatrix>
523  {
524  RCP<ematrix_t> M;
525 
526  try{
527  M = euinput->getUIEpetraCrsMatrix();
528  }
529  catch(std::exception &e){
530  aok = false;
531  std::cout << e.what() << std::endl;
532  }
533  TEST_FAIL_AND_EXIT(*comm, aok, "getEpetraCrsMatrix ", 1);
534 
535  if (rank== 0)
536  std::cout << "Original Epetra matrix" << std::endl;
537 
538  M->Print(std::cout);
539 
540  RCP<const emap_t> emap = Teuchos::rcpFromRef(M->RowMap());
541  RCP<const xemap_t> xmap(new xemap_t(emap));
542 
543  ArrayRCP<zgno_t> newRowIds = roundRobinMap(xmap);
544 
545  zgno_t localNumRows = newRowIds.size();
546 
547  RCP<const ematrix_t> newM;
548  try{
550  localNumRows, newRowIds.getRawPtr());
551  }
552  catch(std::exception &e){
553  aok = false;
554  std::cout << e.what() << std::endl;
555  }
556  TEST_FAIL_AND_EXIT(*comm, aok,
557  " Zoltan2::XpetraTraits<ematrix_t>::doMigration ", 1);
558 
559  if (rank== 0)
560  std::cout << "Migrated Epetra matrix" << std::endl;
561 
562  newM->Print(std::cout);
563  }
564 
565  // XpetraTraits<Epetra_CrsGraph>
566  {
567  RCP<egraph_t> G;
568 
569  try{
570  G = euinput->getUIEpetraCrsGraph();
571  }
572  catch(std::exception &e){
573  aok = false;
574  std::cout << e.what() << std::endl;
575  }
576  TEST_FAIL_AND_EXIT(*comm, aok, "getEpetraCrsGraph ", 1);
577 
578  if (rank== 0)
579  std::cout << "Original Epetra graph" << std::endl;
580 
581  G->Print(std::cout);
582 
583  RCP<const emap_t> emap = Teuchos::rcpFromRef(G->RowMap());
584  RCP<const xemap_t> xmap(new xemap_t(emap));
585  ArrayRCP<zgno_t> newRowIds = roundRobinMap(xmap);
586 
587  zgno_t localNumRows = newRowIds.size();
588 
589  RCP<const egraph_t> newG;
590  try{
592  localNumRows, newRowIds.getRawPtr());
593  }
594  catch(std::exception &e){
595  aok = false;
596  std::cout << e.what() << std::endl;
597  }
598  TEST_FAIL_AND_EXIT(*comm, aok,
599  " Zoltan2::XpetraTraits<egraph_t>::doMigration ", 1);
600 
601  if (rank== 0)
602  std::cout << "Migrated Epetra graph" << std::endl;
603 
604  newG->Print(std::cout);
605  }
606 
607  // XpetraTraits<Epetra_Vector>
608  {
609  RCP<evector_t> V;
610 
611  try{
612  V = rcp(new Epetra_Vector(euinput->getUIEpetraCrsGraph()->RowMap()));
613  V->Random();
614  }
615  catch(std::exception &e){
616  aok = false;
617  std::cout << e.what() << std::endl;
618  }
619  TEST_FAIL_AND_EXIT(*comm, aok, "getEpetraVector", 1);
620 
621  if (rank== 0)
622  std::cout << "Original Epetra vector" << std::endl;
623 
624  V->Print(std::cout);
625 
626  RCP<const emap_t> emap = Teuchos::rcpFromRef(V->Map());
627  RCP<const xemap_t> xmap(new xemap_t(emap));
628  ArrayRCP<zgno_t> newRowIds = roundRobinMap(xmap);
629 
630  zgno_t localNumRows = newRowIds.size();
631 
632  RCP<const evector_t> newV;
633  try{
635  localNumRows, newRowIds.getRawPtr());
636  }
637  catch(std::exception &e){
638  aok = false;
639  std::cout << e.what() << std::endl;
640  }
641  TEST_FAIL_AND_EXIT(*comm, aok,
642  " Zoltan2::XpetraTraits<evector_t>::doMigration ", 1);
643 
644  if (rank== 0)
645  std::cout << "Migrated Epetra vector" << std::endl;
646 
647  newV->Print(std::cout);
648  }
649 
650  // XpetraTraits<Epetra_MultiVector>
651  {
652  RCP<emvector_t> MV;
653 
654  try{
655  MV =
656  rcp(new Epetra_MultiVector(euinput->getUIEpetraCrsGraph()->RowMap(),3));
657  MV->Random();
658  }
659  catch(std::exception &e){
660  aok = false;
661  std::cout << e.what() << std::endl;
662  }
663  TEST_FAIL_AND_EXIT(*comm, aok, "getEpetraMultiVector", 1);
664 
665  if (rank== 0)
666  std::cout << "Original Epetra multivector" << std::endl;
667 
668  MV->Print(std::cout);
669 
670  RCP<const emap_t> emap = Teuchos::rcpFromRef(MV->Map());
671  RCP<const xemap_t> xmap(new xemap_t(emap));
672  ArrayRCP<zgno_t> newRowIds = roundRobinMap(xmap);
673 
674  zgno_t localNumRows = newRowIds.size();
675 
676  RCP<const emvector_t> newMV;
677  try{
679  localNumRows, newRowIds.getRawPtr());
680  }
681  catch(std::exception &e){
682  aok = false;
683  std::cout << e.what() << std::endl;
684  }
685  TEST_FAIL_AND_EXIT(*comm, aok,
686  " Zoltan2::XpetraTraits<emvector_t>::doMigration ", 1);
687 
688  if (rank== 0)
689  std::cout << "Migrated Epetra multivector" << std::endl;
690 
691  newMV->Print(std::cout);
692  }
693 #endif // have epetra data types (int, int, double)
694 
696  // DONE
698 
699  if (rank==0)
700  std::cout << "PASS" << std::endl;
701 }
702 
Tpetra::CrsMatrix< zscalar_t, zlno_t, zgno_t, znode_t > tmatrix_t
Xpetra::CrsMatrix< zscalar_t, zlno_t, zgno_t, znode_t > xmatrix_t
#define TEST_FAIL_AND_EXIT(comm, ok, s, code)
static RCP< User > doMigration(const User &from, size_t numLocalRows, const gno_t *myNewRows)
Migrate the object Given a user object and a new row distribution, create and return a new user objec...
Tpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > tvector_t
Xpetra::CrsGraph< zlno_t, zgno_t, znode_t > xgraph_t
Traits of Xpetra classes, including migration method.
common code used by tests
int main(int argc, char *argv[])
static RCP< User > convertToXpetra(const RCP< User > &a)
Convert the object to its Xpetra wrapped version.
ArrayRCP< zgno_t > roundRobinMap(const RCP< const Xpetra::Map< zlno_t, zgno_t, znode_t > > &m)
int zgno_t
Xpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > xvector_t
Tpetra::CrsGraph< zlno_t, zgno_t, znode_t > tgraph_t
std::string testDataFilePath(".")