Zoltan2 Version of the Day
XpetraTraits.cpp
Go to the documentation of this file.
00001 // @HEADER
00002 //
00003 // ***********************************************************************
00004 //
00005 //   Zoltan2: A package of combinatorial algorithms for scientific computing
00006 //                  Copyright 2012 Sandia Corporation
00007 //
00008 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00009 // the U.S. Government retains certain rights in this software.
00010 //
00011 // Redistribution and use in source and binary forms, with or without
00012 // modification, are permitted provided that the following conditions are
00013 // met:
00014 //
00015 // 1. Redistributions of source code must retain the above copyright
00016 // notice, this list of conditions and the following disclaimer.
00017 //
00018 // 2. Redistributions in binary form must reproduce the above copyright
00019 // notice, this list of conditions and the following disclaimer in the
00020 // documentation and/or other materials provided with the distribution.
00021 //
00022 // 3. Neither the name of the Corporation nor the names of the
00023 // contributors may be used to endorse or promote products derived from
00024 // this software without specific prior written permission.
00025 //
00026 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00027 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00028 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00029 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00030 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00031 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00032 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00033 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00034 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00035 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00036 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00037 //
00038 // Questions? Contact Karen Devine      (kddevin@sandia.gov)
00039 //                    Erik Boman        (egboman@sandia.gov)
00040 //                    Siva Rajamanickam (srajama@sandia.gov)
00041 //
00042 // ***********************************************************************
00043 //
00044 // @HEADER
00045 //
00046 // Basic test of the XpetraTraits definitions.
00047 //
00048 // TODO - a real test would figure out if the migrated objects are
00049 // the same as the original, here we just look at them on stdout.
00050 // TODO look at number of diagonals and max number of entries in
00051 //   Tpetra and Xpetra migrated graphs.  They're garbage.
00052 
00053 #include <Zoltan2_XpetraTraits.hpp>
00054 #include <Zoltan2_TestHelpers.hpp>
00055 
00056 #include <string>
00057 #include <Teuchos_GlobalMPISession.hpp>
00058 #include <Teuchos_DefaultComm.hpp>
00059 #include <Teuchos_RCP.hpp>
00060 #include <Teuchos_Array.hpp>
00061 #include <Teuchos_ArrayRCP.hpp>
00062 #include <Teuchos_Comm.hpp>
00063 #include <Teuchos_VerboseObject.hpp>
00064 #include <Tpetra_CrsMatrix.hpp>
00065 #include <Tpetra_Vector.hpp>
00066 
00067 #include <Xpetra_EpetraUtils.hpp>
00068 #ifdef HAVE_ZOLTAN2_MPI
00069 #include <Epetra_MpiComm.h>
00070 #else
00071 #include <Epetra_SerialComm.h>
00072 #endif
00073 
00074 using namespace std;
00075 using std::string;
00076 using Teuchos::RCP;
00077 using Teuchos::ArrayRCP;
00078 using Teuchos::ArrayView;
00079 using Teuchos::Array;
00080 using Teuchos::rcp;
00081 using Teuchos::Comm;
00082 
00083 ArrayRCP<gno_t> roundRobinMap(
00084     const RCP<const Xpetra::Map<lno_t, gno_t, node_t> > &m)
00085 {
00086   const RCP<const Comm<int> > &comm = m->getComm();
00087   int proc = comm->getRank();
00088   int nprocs = comm->getSize();
00089   gno_t base = m->getMinAllGlobalIndex();
00090   gno_t max = m->getMaxAllGlobalIndex();
00091   size_t globalrows = m->getGlobalNumElements();
00092   if (globalrows != size_t(max - base + 1)){
00093     TEST_FAIL_AND_EXIT(*comm, 0, 
00094       string("Map is invalid for test - fix test"), 1);
00095   }
00096   RCP<Array<gno_t> > mygids = rcp(new Array<gno_t>);
00097   gno_t firstgno_t = proc; 
00098   if (firstgno_t < base){
00099     gno_t n = base % proc;
00100     if (n>0)
00101       firstgno_t = base - n + proc;
00102     else
00103       firstgno_t = base;
00104   }
00105   for (gno_t gid=firstgno_t; gid <= max; gid+=nprocs){
00106     (*mygids).append(gid);
00107   }
00108 
00109   ArrayRCP<gno_t> newIdArcp = Teuchos::arcp(mygids);
00110 
00111   return newIdArcp;
00112 }
00113 
00114 int main(int argc, char *argv[])
00115 {
00116   Teuchos::GlobalMPISession session(&argc, &argv);
00117   RCP<const Comm<int> > comm = Teuchos::DefaultComm<int>::getComm();
00118   int rank = comm->getRank();
00119 
00120   Teuchos::RCP<Teuchos::FancyOStream> outStream = 
00121     Teuchos::VerboseObjectBase::getDefaultOStream();
00122   Teuchos::EVerbosityLevel v=Teuchos::VERB_EXTREME;
00123 
00124   typedef UserInputForTests uinput_t;
00125   typedef Tpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t> tmatrix_t;
00126   typedef Tpetra::CrsGraph<lno_t,gno_t,node_t> tgraph_t;
00127   typedef Tpetra::Vector<scalar_t,lno_t,gno_t,node_t> tvector_t;
00128   typedef Tpetra::MultiVector<scalar_t,lno_t,gno_t,node_t> tmvector_t;
00129   typedef Xpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t> xmatrix_t;
00130   typedef Xpetra::CrsGraph<lno_t,gno_t,node_t> xgraph_t;
00131   typedef Xpetra::Vector<scalar_t,lno_t,gno_t,node_t> xvector_t;
00132   typedef Xpetra::MultiVector<scalar_t,lno_t,gno_t,node_t> xmvector_t;
00133   typedef Xpetra::TpetraMap<lno_t,gno_t,node_t> xtmap_t;
00134 
00135   // Create object that can give us test Tpetra and Xpetra input.
00136 
00137   RCP<uinput_t> uinput;
00138 
00139   try{
00140     uinput = 
00141       rcp(new uinput_t(testDataFilePath,std::string("simple"), comm, true));
00142   }
00143   catch(std::exception &e){
00144     TEST_FAIL_AND_EXIT(*comm, 0, string("input ")+e.what(), 1);
00145   }
00146 
00148   //   Tpetra::CrsMatrix
00149   //   Tpetra::CrsGraph
00150   //   Tpetra::Vector
00151   //   Tpetra::MultiVector
00153 
00154   // XpetraTraits<Tpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> > 
00155   {
00156     RCP<tmatrix_t> M;
00157   
00158     try{
00159       M = uinput->getTpetraCrsMatrix();
00160     }
00161     catch(std::exception &e){
00162       TEST_FAIL_AND_EXIT(*comm, 0, 
00163         string("getTpetraCrsMatrix ")+e.what(), 1);
00164     }
00165   
00166     if (rank== 0)
00167       std::cout << "Original Tpetra matrix " << M->getGlobalNumRows()
00168         << " x " << M->getGlobalNumCols() << std::endl;
00169 
00170     M->describe(*outStream,v);
00171 
00172     RCP<const xtmap_t> xmap(new xtmap_t(M->getRowMap()));
00173 
00174     ArrayRCP<gno_t> newRowIds = roundRobinMap(xmap);
00175   
00176     gno_t localNumRows = newRowIds.size();
00177   
00178     RCP<const tmatrix_t> newM;
00179     try{
00180       newM = Zoltan2::XpetraTraits<tmatrix_t>::doMigration(
00181         rcp_const_cast<const tmatrix_t>(M),
00182         localNumRows, newRowIds.getRawPtr());
00183     }
00184     catch(std::exception &e){
00185       TEST_FAIL_AND_EXIT(*comm, 0, 
00186         string(" Zoltan2::XpetraTraits<tmatrix_t>::doMigration ")+e.what(), 1);
00187     }
00188 
00189     if (rank== 0)
00190       std::cout << "Migrated Tpetra matrix" << std::endl;
00191   
00192     newM->describe(*outStream,v);
00193   }
00194 
00195   // XpetraTraits<Tpetra::CrsGraph<scalar_t, lno_t, gno_t, node_t> > 
00196   {
00197     RCP<tgraph_t> G;
00198   
00199     try{
00200       G = uinput->getTpetraCrsGraph();
00201     }
00202     catch(std::exception &e){
00203       TEST_FAIL_AND_EXIT(*comm, 0, 
00204         string("getTpetraCrsGraph ")+e.what(), 1);
00205     }
00206   
00207     if (rank== 0)
00208       std::cout << "Original Tpetra graph" << std::endl;
00209   
00210     G->describe(*outStream,v);
00211   
00212     RCP<const xtmap_t> xmap(new xtmap_t(G->getRowMap()));
00213     ArrayRCP<gno_t> newRowIds = roundRobinMap(xmap);
00214   
00215     gno_t localNumRows = newRowIds.size();
00216   
00217     RCP<const tgraph_t> newG;
00218     try{
00219       newG = Zoltan2::XpetraTraits<tgraph_t>::doMigration(
00220         rcp_const_cast<const tgraph_t>(G),
00221         localNumRows, newRowIds.getRawPtr());
00222     }
00223     catch(std::exception &e){
00224       TEST_FAIL_AND_EXIT(*comm, 0, 
00225         string(" Zoltan2::XpetraTraits<tgraph_t>::doMigration ")+e.what(), 1);
00226     }
00227   
00228     if (rank== 0)
00229       std::cout << "Migrated Tpetra graph" << std::endl;
00230   
00231     newG->describe(*outStream,v);
00232   }
00233 
00234   // XpetraTraits<Tpetra::Vector<scalar_t, lno_t, gno_t, node_t>> 
00235   {
00236     RCP<tvector_t> V;
00237   
00238     try{
00239       V = uinput->getTpetraVector();
00240     }
00241     catch(std::exception &e){
00242       TEST_FAIL_AND_EXIT(*comm, 0, 
00243         string("getTpetraVector")+e.what(), 1);
00244     }
00245   
00246     if (rank== 0)
00247       std::cout << "Original Tpetra vector" << std::endl;
00248   
00249     V->describe(*outStream,v);
00250   
00251     RCP<const xtmap_t> xmap(new xtmap_t(V->getMap()));
00252     ArrayRCP<gno_t> newRowIds = roundRobinMap(xmap);
00253   
00254     gno_t localNumRows = newRowIds.size();
00255   
00256     RCP<const tvector_t> newV;
00257     try{
00258       newV = Zoltan2::XpetraTraits<tvector_t>::doMigration(
00259         rcp_const_cast<const tvector_t>(V),
00260         localNumRows, newRowIds.getRawPtr());
00261     }
00262     catch(std::exception &e){
00263       TEST_FAIL_AND_EXIT(*comm, 0, 
00264         string(" Zoltan2::XpetraTraits<tvector_t>::doMigration ")+e.what(), 1);
00265     }
00266   
00267     if (rank== 0)
00268       std::cout << "Migrated Tpetra vector" << std::endl;
00269   
00270     newV->describe(*outStream,v);
00271   }
00272 
00273   // XpetraTraits<Tpetra::MultiVector<scalar_t, lno_t, gno_t, node_t>> 
00274   {
00275     RCP<tmvector_t> MV;
00276   
00277     try{
00278       MV = uinput->getTpetraMultiVector(3);
00279     }
00280     catch(std::exception &e){
00281       TEST_FAIL_AND_EXIT(*comm, 0, 
00282         string("getTpetraMultiVector")+e.what(), 1);
00283     }
00284   
00285     if (rank== 0)
00286       std::cout << "Original Tpetra multivector" << std::endl;
00287   
00288     MV->describe(*outStream,v);
00289   
00290     RCP<const xtmap_t> xmap(new xtmap_t(MV->getMap()));
00291     ArrayRCP<gno_t> newRowIds = roundRobinMap(xmap);
00292   
00293     gno_t localNumRows = newRowIds.size();
00294   
00295     RCP<const tmvector_t> newMV;
00296     try{
00297       newMV = Zoltan2::XpetraTraits<tmvector_t>::doMigration(
00298         rcp_const_cast<const tmvector_t>(MV),
00299         localNumRows, newRowIds.getRawPtr());
00300     }
00301     catch(std::exception &e){
00302       TEST_FAIL_AND_EXIT(*comm, 0, 
00303         string(" Zoltan2::XpetraTraits<tmvector_t>::doMigration ")+e.what(), 1);
00304     }
00305   
00306     if (rank== 0)
00307       std::cout << "Migrated Tpetra multivector" << std::endl;
00308   
00309     newMV->describe(*outStream,v);
00310   }
00311 
00313   //   Xpetra::CrsMatrix
00314   //   Xpetra::CrsGraph
00315   //   Xpetra::Vector
00316   //   Xpetra::MultiVector
00318 
00319   // XpetraTraits<Xpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> > 
00320   {
00321     RCP<xmatrix_t> M;
00322   
00323     try{
00324       M = uinput->getXpetraCrsMatrix();
00325     }
00326     catch(std::exception &e){
00327       TEST_FAIL_AND_EXIT(*comm, 0, 
00328         string("getXpetraCrsMatrix ")+e.what(), 1);
00329     }
00330   
00331     if (rank== 0)
00332       std::cout << "Original Xpetra matrix" << std::endl;
00333   
00334     M->describe(*outStream,v);
00335   
00336     ArrayRCP<gno_t> newRowIds = roundRobinMap(M->getRowMap());
00337   
00338     gno_t localNumRows = newRowIds.size();
00339   
00340     RCP<const xmatrix_t> newM;
00341     try{
00342       newM = Zoltan2::XpetraTraits<xmatrix_t>::doMigration(
00343         rcp_const_cast<const xmatrix_t>(M),
00344         localNumRows, newRowIds.getRawPtr());
00345     }
00346     catch(std::exception &e){
00347       TEST_FAIL_AND_EXIT(*comm, 0, 
00348         string(" Zoltan2::XpetraTraits<xmatrix_t>::doMigration ")+e.what(), 1);
00349     }
00350 
00351     if (rank== 0)
00352       std::cout << "Migrated Xpetra matrix" << std::endl;
00353   
00354     newM->describe(*outStream,v);
00355   }
00356 
00357   // XpetraTraits<Xpetra::CrsGraph<scalar_t, lno_t, gno_t, node_t> > 
00358   {
00359     RCP<xgraph_t> G;
00360   
00361     try{
00362       G = uinput->getXpetraCrsGraph();
00363     }
00364     catch(std::exception &e){
00365       TEST_FAIL_AND_EXIT(*comm, 0, 
00366         string("getXpetraCrsGraph ")+e.what(), 1);
00367     }
00368   
00369     if (rank== 0)
00370       std::cout << "Original Xpetra graph" << std::endl;
00371   
00372     G->describe(*outStream,v);
00373   
00374     ArrayRCP<gno_t> newRowIds = roundRobinMap(G->getRowMap());
00375   
00376     gno_t localNumRows = newRowIds.size();
00377   
00378     RCP<const xgraph_t> newG;
00379     try{
00380       newG = Zoltan2::XpetraTraits<xgraph_t>::doMigration(
00381         rcp_const_cast<const xgraph_t>(G),
00382         localNumRows, newRowIds.getRawPtr());
00383     }
00384     catch(std::exception &e){
00385       TEST_FAIL_AND_EXIT(*comm, 0, 
00386         string(" Zoltan2::XpetraTraits<xgraph_t>::doMigration ")+e.what(), 1);
00387     }
00388   
00389     if (rank== 0)
00390       std::cout << "Migrated Xpetra graph" << std::endl;
00391   
00392     newG->describe(*outStream,v);
00393   }
00394 
00395   // XpetraTraits<Xpetra::Vector<scalar_t, lno_t, gno_t, node_t>> 
00396   {
00397     RCP<xvector_t> V;
00398   
00399     try{
00400       V = uinput->getXpetraVector();
00401     }
00402     catch(std::exception &e){
00403       TEST_FAIL_AND_EXIT(*comm, 0, 
00404         string("getXpetraVector")+e.what(), 1);
00405     }
00406   
00407     if (rank== 0)
00408       std::cout << "Original Xpetra vector" << std::endl;
00409   
00410     V->describe(*outStream,v);
00411   
00412     ArrayRCP<gno_t> newRowIds = roundRobinMap(V->getMap());
00413   
00414     gno_t localNumRows = newRowIds.size();
00415   
00416     RCP<const xvector_t> newV;
00417     try{
00418       newV = Zoltan2::XpetraTraits<xvector_t>::doMigration(
00419         rcp_const_cast<const xvector_t>(V),
00420         localNumRows, newRowIds.getRawPtr());
00421     }
00422     catch(std::exception &e){
00423       TEST_FAIL_AND_EXIT(*comm, 0, 
00424         string(" Zoltan2::XpetraTraits<xvector_t>::doMigration ")+e.what(), 1);
00425     }
00426   
00427     if (rank== 0)
00428       std::cout << "Migrated Xpetra vector" << std::endl;
00429   
00430     newV->describe(*outStream,v);
00431   }
00432 
00433   // XpetraTraits<Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t>> 
00434   {
00435     RCP<xmvector_t> MV;
00436   
00437     try{
00438       MV = uinput->getXpetraMultiVector(3);
00439     }
00440     catch(std::exception &e){
00441       TEST_FAIL_AND_EXIT(*comm, 0, 
00442         string("getXpetraMultiVector")+e.what(), 1);
00443     }
00444   
00445     if (rank== 0)
00446       std::cout << "Original Xpetra multivector" << std::endl;
00447   
00448     MV->describe(*outStream,v);
00449   
00450     ArrayRCP<gno_t> newRowIds = roundRobinMap(MV->getMap());
00451   
00452     gno_t localNumRows = newRowIds.size();
00453   
00454     RCP<const xmvector_t> newMV;
00455     try{
00456       newMV = Zoltan2::XpetraTraits<xmvector_t>::doMigration(
00457         rcp_const_cast<const xmvector_t>(MV),
00458         localNumRows, newRowIds.getRawPtr());
00459     }
00460     catch(std::exception &e){
00461       TEST_FAIL_AND_EXIT(*comm, 0, 
00462         string(" Zoltan2::XpetraTraits<xmvector_t>::doMigration ")+e.what(), 1);
00463     }
00464   
00465     if (rank== 0)
00466       std::cout << "Migrated Xpetra multivector" << std::endl;
00467   
00468     newMV->describe(*outStream,v);
00469   }
00470 
00471 #ifdef HAVE_EPETRA_DATA_TYPES
00472 
00473   //   Epetra_CrsMatrix
00474   //   Epetra_CrsGraph
00475   //   Epetra_Vector
00476   //   Epetra_MultiVector
00478 
00479   typedef Epetra_CrsMatrix ematrix_t;
00480   typedef Epetra_CrsGraph egraph_t;
00481   typedef Epetra_Vector evector_t;
00482   typedef Epetra_MultiVector emvector_t;
00483   typedef Xpetra::EpetraMap xemap_t;
00484   typedef Epetra_BlockMap emap_t;
00485 
00486   // Create object that can give us test Epetra input.
00487 
00488   RCP<uinput_t> euinput;
00489 
00490   try{
00491     euinput = 
00492       rcp(new uinput_t(testDataFilePath,std::string("simple"), comm, true));
00493   }
00494   catch(std::exception &e){
00495     TEST_FAIL_AND_EXIT(*comm, 0, string("epetra input ")+e.what(), 1);
00496   }
00497 
00498   // XpetraTraits<Epetra_CrsMatrix> 
00499   {
00500     RCP<ematrix_t> M;
00501   
00502     try{
00503       M = euinput->getEpetraCrsMatrix();
00504     }
00505     catch(std::exception &e){
00506       TEST_FAIL_AND_EXIT(*comm, 0, 
00507         string("getEpetraCrsMatrix ")+e.what(), 1);
00508     }
00509 
00510     if (rank== 0)
00511       std::cout << "Original Epetra matrix" << std::endl;
00512   
00513     M->Print(std::cout);
00514   
00515     RCP<const emap_t> emap = Teuchos::rcpFromRef(M->RowMap());
00516     RCP<const xemap_t> xmap(new xemap_t(emap));
00517 
00518     ArrayRCP<gno_t> newRowIds = roundRobinMap(xmap);
00519   
00520     gno_t localNumRows = newRowIds.size();
00521   
00522     RCP<const ematrix_t> newM;
00523     try{
00524       newM = Zoltan2::XpetraTraits<ematrix_t>::doMigration(
00525         rcp_const_cast<const ematrix_t>(M),
00526         localNumRows, newRowIds.getRawPtr());
00527     }
00528     catch(std::exception &e){
00529       TEST_FAIL_AND_EXIT(*comm, 0, 
00530         string(" Zoltan2::XpetraTraits<ematrix_t>::doMigration ")+e.what(), 1);
00531     }
00532 
00533     if (rank== 0)
00534       std::cout << "Migrated Epetra matrix" << std::endl;
00535   
00536     newM->Print(std::cout);
00537   }
00538 
00539   // XpetraTraits<Epetra_CrsGraph> 
00540   {
00541     RCP<egraph_t> G;
00542   
00543     try{
00544       G = euinput->getEpetraCrsGraph();
00545     }
00546     catch(std::exception &e){
00547       TEST_FAIL_AND_EXIT(*comm, 0, 
00548         string("getEpetraCrsGraph ")+e.what(), 1);
00549     }
00550   
00551     if (rank== 0)
00552       std::cout << "Original Epetra graph" << std::endl;
00553   
00554     G->Print(std::cout);
00555   
00556     RCP<const emap_t> emap = Teuchos::rcpFromRef(G->RowMap());
00557     RCP<const xemap_t> xmap(new xemap_t(emap));
00558     ArrayRCP<gno_t> newRowIds = roundRobinMap(xmap);
00559   
00560     gno_t localNumRows = newRowIds.size();
00561   
00562     RCP<const egraph_t> newG;
00563     try{
00564       newG = Zoltan2::XpetraTraits<egraph_t>::doMigration(
00565         rcp_const_cast<const egraph_t>(G),
00566         localNumRows, newRowIds.getRawPtr());
00567     }
00568     catch(std::exception &e){
00569       TEST_FAIL_AND_EXIT(*comm, 0, 
00570         string(" Zoltan2::XpetraTraits<egraph_t>::doMigration ")+e.what(), 1);
00571     }
00572   
00573     if (rank== 0)
00574       std::cout << "Migrated Epetra graph" << std::endl;
00575   
00576     newG->Print(std::cout);
00577   }
00578 
00579   // XpetraTraits<Epetra_Vector>
00580   {
00581     RCP<evector_t> V;
00582   
00583     try{
00584       V = euinput->getEpetraVector();
00585     }
00586     catch(std::exception &e){
00587       TEST_FAIL_AND_EXIT(*comm, 0, 
00588         string("getEpetraVector")+e.what(), 1);
00589     }
00590   
00591     if (rank== 0)
00592       std::cout << "Original Epetra vector" << std::endl;
00593   
00594     V->Print(std::cout);
00595   
00596     RCP<const emap_t> emap = Teuchos::rcpFromRef(V->Map());
00597     RCP<const xemap_t> xmap(new xemap_t(emap));
00598     ArrayRCP<gno_t> newRowIds = roundRobinMap(xmap);
00599   
00600     gno_t localNumRows = newRowIds.size();
00601   
00602     RCP<const evector_t> newV;
00603     try{
00604       newV = Zoltan2::XpetraTraits<evector_t>::doMigration(
00605         rcp_const_cast<const evector_t>(V),
00606         localNumRows, newRowIds.getRawPtr());
00607     }
00608     catch(std::exception &e){
00609       TEST_FAIL_AND_EXIT(*comm, 0, 
00610         string(" Zoltan2::XpetraTraits<evector_t>::doMigration ")+e.what(), 1);
00611     }
00612   
00613     if (rank== 0)
00614       std::cout << "Migrated Epetra vector" << std::endl;
00615   
00616     newV->Print(std::cout);
00617   }
00618 
00619   // XpetraTraits<Epetra_MultiVector>
00620   {
00621     RCP<emvector_t> MV;
00622   
00623     try{
00624       MV = euinput->getEpetraMultiVector(3);
00625     }
00626     catch(std::exception &e){
00627       TEST_FAIL_AND_EXIT(*comm, 0, 
00628         string("getEpetraMultiVector")+e.what(), 1);
00629     }
00630   
00631     if (rank== 0)
00632       std::cout << "Original Epetra multivector" << std::endl;
00633   
00634     MV->Print(std::cout);
00635   
00636     RCP<const emap_t> emap = Teuchos::rcpFromRef(MV->Map());
00637     RCP<const xemap_t> xmap(new xemap_t(emap));
00638     ArrayRCP<gno_t> newRowIds = roundRobinMap(xmap);
00639   
00640     gno_t localNumRows = newRowIds.size();
00641   
00642     RCP<const emvector_t> newMV;
00643     try{
00644       newMV = Zoltan2::XpetraTraits<emvector_t>::doMigration(
00645         rcp_const_cast<const emvector_t>(MV),
00646         localNumRows, newRowIds.getRawPtr());
00647     }
00648     catch(std::exception &e){
00649       TEST_FAIL_AND_EXIT(*comm, 0, 
00650         string(" Zoltan2::XpetraTraits<emvector_t>::doMigration ")+e.what(), 1);
00651     }
00652   
00653     if (rank== 0)
00654       std::cout << "Migrated Epetra multivector" << std::endl;
00655   
00656     newMV->Print(std::cout);
00657   }
00658 #endif   // have epetra data types (int, int, double)
00659 
00661   // DONE
00663 
00664   if (rank==0)
00665     std::cout << "PASS" << std::endl;
00666 }
00667