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