Zoltan2 Version of the Day
Zoltan2_XpetraTraits.hpp
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 
00050 #ifndef _ZOLTAN2_XPETRATRAITS_HPP_
00051 #define _ZOLTAN2_XPETRATRAITS_HPP_
00052 
00053 #include <Zoltan2_InputTraits.hpp>
00054 #include <Zoltan2_Standards.hpp>
00055 
00056 #include <Xpetra_EpetraCrsMatrix.hpp>
00057 #include <Xpetra_TpetraCrsMatrix.hpp>
00058 #include <Xpetra_TpetraRowMatrix.hpp>
00059 #include <Xpetra_EpetraVector.hpp>
00060 #include <Xpetra_TpetraVector.hpp>
00061 #include <Xpetra_EpetraUtils.hpp>
00062 #include <Tpetra_Vector.hpp>
00063 
00064 namespace Zoltan2 {
00065 
00067 // Extra traits needed only for Xpetra matrices and graphs
00068 
00091 template <typename User>
00092 struct XpetraTraits
00093 {
00096   static inline RCP<const User> convertToXpetra(const RCP<const User> &a) 
00097   {
00098     return a;
00099   }
00100 
00103   typedef long gno_t;
00104 
00107   typedef int lno_t;
00108 
00114   static RCP<const User> doMigration(const RCP<const User> &from,
00115       lno_t numLocalRows, const gno_t *myNewRows)
00116   {
00117     return from;
00118   }
00119 };
00120 
00121 #ifndef DOXYGEN_SHOULD_SKIP_THIS
00122 
00124 // Tpetra::CrsMatrix
00125 template <typename scalar_t,
00126           typename lno_t,
00127           typename gno_t,
00128           typename node_t>
00129 struct XpetraTraits<Tpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> >
00130 {
00131   typedef typename Xpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t> 
00132                    xmatrix_t;
00133   typedef typename Xpetra::TpetraCrsMatrix<scalar_t,lno_t,gno_t,node_t> 
00134                    xtmatrix_t;
00135   typedef typename Tpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t> 
00136                    tmatrix_t;
00137 
00138   static inline RCP<const xmatrix_t> convertToXpetra(
00139     const RCP<const tmatrix_t> &a)
00140     {
00141       return rcp(new xtmatrix_t(rcp_const_cast<tmatrix_t>(a)));
00142     }
00143 
00144   static RCP<const tmatrix_t> doMigration(const RCP<const tmatrix_t> &from,
00145       lno_t numLocalRows, const gno_t *myNewRows)
00146   {
00147     typedef Tpetra::Map<lno_t, gno_t, node_t> map_t;
00148     lno_t base = 0;
00149 
00150     // source map
00151     const RCP<const map_t> &smap = from->getRowMap();
00152     int oldNumElts = smap->getNodeNumElements();
00153     gno_t numGlobalRows = smap->getGlobalNumElements();
00154 
00155     // target map
00156     ArrayView<const gno_t> rowList(myNewRows, numLocalRows);
00157     const RCP<const Teuchos::Comm<int> > &comm = from->getComm();
00158     RCP<const map_t> tmap = rcp(
00159       new map_t(numGlobalRows, rowList, base, comm));
00160     int newNumElts = numLocalRows;
00161 
00162     // importer
00163     Tpetra::Import<lno_t, gno_t, node_t> importer(smap, tmap);
00164 
00165     // number of non zeros in my new rows
00166     typedef Tpetra::Vector<scalar_t, lno_t, gno_t, node_t> vector_t;
00167     vector_t numOld(smap);  // TODO These vectors should have scalar = size_t, 
00168     vector_t numNew(tmap);  // but explicit instantiation does not yet support that.
00169     for (int lid=0; lid < oldNumElts; lid++){
00170       numOld.replaceGlobalValue(smap->getGlobalElement(lid), 
00171         scalar_t(from->getNumEntriesInLocalRow(lid)));
00172     }
00173     numNew.doImport(numOld, importer, Tpetra::INSERT);
00174 
00175     // TODO Could skip this copy if could declare vector with scalar = size_t.
00176     ArrayRCP<size_t> nnz(newNumElts);
00177     if (newNumElts > 0){
00178       ArrayRCP<scalar_t> ptr = numNew.getDataNonConst(0);
00179       for (int lid=0; lid < newNumElts; lid++){
00180         nnz[lid] = static_cast<size_t>(ptr[lid]);
00181       }
00182     }
00183 
00184     // target matrix
00185     RCP<tmatrix_t> M = rcp(new tmatrix_t(tmap, nnz, Tpetra::StaticProfile));
00186     M->doImport(*from, importer, Tpetra::INSERT);
00187     M->fillComplete();
00188 
00189     return rcp_const_cast<const tmatrix_t>(M);
00190   }
00191 };
00192 
00194 // Epetra_CrsMatrix
00195 
00196 template <>
00197 struct XpetraTraits<Epetra_CrsMatrix>
00198 {
00199   typedef InputTraits<Epetra_CrsMatrix>::scalar_t scalar_t;
00200   typedef InputTraits<Epetra_CrsMatrix>::lno_t lno_t;
00201   typedef InputTraits<Epetra_CrsMatrix>::gno_t gno_t;
00202   typedef InputTraits<Epetra_CrsMatrix>::node_t node_t;
00203 
00204   static inline RCP<const Xpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t> >
00205     convertToXpetra(const RCP<const Epetra_CrsMatrix> &a)
00206     {
00207       return rcp(new Xpetra::EpetraCrsMatrix(
00208                              rcp_const_cast<Epetra_CrsMatrix>(a)));
00209     }
00210 
00211 
00212   static RCP<Epetra_CrsMatrix> doMigration(
00213       const RCP<const Epetra_CrsMatrix> &from,
00214       lno_t numLocalRows, const gno_t *myNewRows)
00215   {
00216     lno_t base = 0;
00217 
00218     // source map
00219     const Epetra_Map &smap = from->RowMap();
00220     int oldNumElts = smap.NumMyElements();
00221     gno_t numGlobalRows = smap.NumGlobalElements();
00222 
00223     // target map
00224     const Epetra_Comm &comm = from->Comm();
00225     Epetra_Map tmap(numGlobalRows, numLocalRows, myNewRows, base, comm);
00226     int newNumElts = numLocalRows;
00227 
00228     // importer
00229     Epetra_Import importer(tmap, smap);
00230 
00231     // number of non zeros in my new rows
00232     Epetra_Vector numOld(smap);
00233     Epetra_Vector numNew(tmap);
00234 
00235     for (int lid=0; lid < oldNumElts; lid++){
00236       numOld[lid] = from->NumMyEntries(lid);
00237     }
00238     numNew.Import(numOld, importer, Insert);
00239 
00240     Array<int> nnz(newNumElts);
00241     for (int lid=0; lid < newNumElts; lid++){
00242       nnz[lid] = static_cast<int>(numNew[lid]);
00243     }
00244 
00245     // target matrix
00246     RCP<Epetra_CrsMatrix> M = rcp(
00247       new Epetra_CrsMatrix(Copy, tmap, nnz.getRawPtr(), true));
00248     M->Import(*from, importer, Insert);
00249     M->FillComplete();
00250 
00251     return M;
00252   }
00253 };
00254 
00256 // Xpetra::CrsMatrix
00257 // KDDKDD:  Do we need specializations for Xpetra::EpetraCrsMatrix and
00258 // KDDKDD:  Xpetra::TpetraCrsMatrix
00259 template <typename scalar_t,
00260           typename lno_t,
00261           typename gno_t,
00262           typename node_t>
00263 struct XpetraTraits<Xpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> >
00264 {
00265   typedef Xpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> x_matrix_t;
00266   typedef Xpetra::TpetraCrsMatrix<scalar_t, lno_t, gno_t, node_t> xt_matrix_t;
00267   typedef Tpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t> t_matrix_t; 
00268 
00269   static inline RCP<const x_matrix_t>
00270     convertToXpetra( const RCP<const x_matrix_t > &a)
00271     {
00272       return a;
00273     }
00274 
00275   static RCP<const x_matrix_t> doMigration(const RCP<const x_matrix_t> &from,
00276       lno_t numLocalRows, const gno_t *myNewRows)
00277   {
00278     Xpetra::UnderlyingLib lib = from->getRowMap()->lib();
00279 
00280     if (lib == Xpetra::UseEpetra){
00281        throw logic_error("compiler should have used specialization");
00282     } else{
00283       // Do the import with the Tpetra::CrsMatrix traits object
00284       const x_matrix_t *xm = from.get();
00285       const xt_matrix_t *xtm = dynamic_cast<const xt_matrix_t *>(xm);
00286       RCP<const t_matrix_t> tm = xtm->getTpetra_CrsMatrix();
00287 
00288       RCP<const t_matrix_t> tmnew = XpetraTraits<t_matrix_t>::doMigration(
00289         tm, numLocalRows, myNewRows);
00290 
00291       RCP<const x_matrix_t> xmnew = 
00292         XpetraTraits<t_matrix_t>::convertToXpetra(tmnew);
00293 
00294       return xmnew;
00295     }
00296   }
00297 };
00298 
00300 // Xpetra::CrsMatrix specialization
00301 
00302 template <typename node_t>
00303 struct XpetraTraits<Xpetra::CrsMatrix<double, int, int, node_t> >
00304 {
00305   typedef double scalar_t;
00306   typedef int lno_t;
00307   typedef int gno_t;
00308   typedef Xpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> x_matrix_t;
00309   typedef Xpetra::TpetraCrsMatrix<scalar_t, lno_t, gno_t, node_t> xt_matrix_t;
00310   typedef Tpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t> t_matrix_t; 
00311   typedef Xpetra::EpetraCrsMatrix xe_matrix_t;
00312   typedef Epetra_CrsMatrix e_matrix_t; 
00313 
00314   static inline RCP<const x_matrix_t>
00315     convertToXpetra( const RCP<const x_matrix_t > &a)
00316     {
00317       return a;
00318     }
00319 
00320   static RCP<const x_matrix_t> doMigration(const RCP<const x_matrix_t> &from,
00321       lno_t numLocalRows, const gno_t *myNewRows)
00322   {
00323     Xpetra::UnderlyingLib lib = from->getRowMap()->lib();
00324     const x_matrix_t *xm = from.get();
00325 
00326     if (lib == Xpetra::UseEpetra){
00327       // Do the import with the Epetra_CrsMatrix traits object
00328       const xe_matrix_t *xem = dynamic_cast<const xe_matrix_t *>(xm);
00329       RCP<const e_matrix_t> em = xem->getEpetra_CrsMatrix();
00330 
00331       RCP<const e_matrix_t> emnew = XpetraTraits<e_matrix_t>::doMigration(
00332         em, numLocalRows, myNewRows);
00333 
00334       RCP<const x_matrix_t> xmnew = 
00335         XpetraTraits<e_matrix_t>::convertToXpetra(emnew);
00336 
00337       return xmnew;
00338     } else{
00339       // Do the import with the Tpetra::CrsMatrix traits object
00340       const xt_matrix_t *xtm = dynamic_cast<const xt_matrix_t *>(xm);
00341       RCP<const t_matrix_t> tm = xtm->getTpetra_CrsMatrix();
00342 
00343       RCP<const t_matrix_t> tmnew = XpetraTraits<t_matrix_t>::doMigration(
00344         tm, numLocalRows, myNewRows);
00345 
00346       RCP<const x_matrix_t> xmnew = 
00347         XpetraTraits<t_matrix_t>::convertToXpetra(tmnew);
00348 
00349       return xmnew;
00350     }
00351   }
00352 };
00353 
00354 
00356 // Tpetra::CrsGraph
00357 template <typename lno_t,
00358           typename gno_t,
00359           typename node_t>
00360 struct XpetraTraits<Tpetra::CrsGraph<lno_t, gno_t, node_t> >
00361 {
00362   typedef typename Xpetra::CrsGraph<lno_t, gno_t, node_t> xgraph_t;
00363   typedef typename Xpetra::TpetraCrsGraph<lno_t, gno_t, node_t> xtgraph_t;
00364   typedef typename Tpetra::CrsGraph<lno_t, gno_t, node_t> tgraph_t;
00365 
00366   static inline RCP<const xgraph_t> convertToXpetra(
00367     const RCP<const tgraph_t> &a)
00368     {
00369       return rcp(new xtgraph_t(rcp_const_cast<tgraph_t>(a)));
00370     }
00371 
00372   static RCP<const tgraph_t> doMigration(const RCP<const tgraph_t> &from,
00373       lno_t numLocalRows, const gno_t *myNewRows)
00374   {
00375     typedef Tpetra::Map<lno_t, gno_t, node_t> map_t;
00376     lno_t base = 0;
00377 
00378     // source map
00379     const RCP<const map_t> &smap = from->getRowMap();
00380     int oldNumElts = smap->getNodeNumElements();
00381     gno_t numGlobalRows = smap->getGlobalNumElements();
00382 
00383     // target map
00384     ArrayView<const gno_t> rowList(myNewRows, numLocalRows);
00385     const RCP<const Teuchos::Comm<int> > &comm = from->getComm();
00386     RCP<const map_t> tmap = rcp(
00387       new map_t(numGlobalRows, rowList, base, comm));
00388 
00389     // importer
00390     Tpetra::Import<lno_t, gno_t, node_t> importer(smap, tmap);
00391 
00392     // number of entries in my new rows
00393     typedef Tpetra::Vector<gno_t, lno_t, gno_t, node_t> vector_t;
00394     vector_t numOld(smap);
00395     vector_t numNew(tmap);
00396     for (int lid=0; lid < oldNumElts; lid++){
00397       numOld.replaceGlobalValue(smap->getGlobalElement(lid),
00398         from->getNumEntriesInLocalRow(lid));
00399     }
00400     numNew.doImport(numOld, importer, Tpetra::INSERT);
00401 
00402     size_t numElts = tmap->getNodeNumElements();
00403     ArrayRCP<const gno_t> nnz;
00404     if (numElts > 0)
00405       nnz = numNew.getData(0);    // hangs if vector len == 0
00406 
00407     ArrayRCP<const size_t> nnz_size_t;
00408 
00409     if (numElts && sizeof(gno_t) != sizeof(size_t)){
00410       size_t *vals = new size_t [numElts];
00411       nnz_size_t = arcp(vals, 0, numElts, true);
00412       for (size_t i=0; i < numElts; i++){
00413         vals[i] = static_cast<size_t>(nnz[i]);
00414       }
00415     }
00416     else{
00417       nnz_size_t = arcp_reinterpret_cast<const size_t>(nnz);
00418     }
00419 
00420     // target graph
00421     RCP<tgraph_t> G = rcp(new tgraph_t(tmap, nnz_size_t, Tpetra::StaticProfile));
00422 
00423     G->doImport(*from, importer, Tpetra::INSERT);
00424     G->fillComplete();
00425 
00426     return G;
00427   }
00428 
00429 };
00430 
00432 // Epetra_CrsGraph
00433 template < >
00434 struct XpetraTraits<Epetra_CrsGraph>
00435 {
00436   typedef InputTraits<Epetra_CrsGraph>::lno_t    lno_t;
00437   typedef InputTraits<Epetra_CrsGraph>::gno_t    gno_t;
00438   typedef InputTraits<Epetra_CrsGraph>::node_t   node_t;
00439   static inline RCP<const Xpetra::CrsGraph<lno_t,gno_t,node_t> >
00440     convertToXpetra(const RCP<const Epetra_CrsGraph> &a)
00441     {
00442       return rcp(new Xpetra::EpetraCrsGraph(
00443                              rcp_const_cast<Epetra_CrsGraph>(a)));
00444     }
00445 
00446   static RCP<const Epetra_CrsGraph> doMigration(
00447       const RCP<const Epetra_CrsGraph> &from,
00448       lno_t numLocalRows, const gno_t *myNewRows)
00449   {
00450     lno_t base = 0;
00451 
00452     // source map
00453     const Epetra_BlockMap &smap = from->RowMap();
00454     gno_t numGlobalRows = smap.NumGlobalElements();
00455     lno_t oldNumElts = smap.NumMyElements();
00456 
00457     // target map
00458     const Epetra_Comm &comm = from->Comm();
00459     Epetra_BlockMap tmap(numGlobalRows, numLocalRows, 
00460        myNewRows, 1, base, comm);
00461     lno_t newNumElts = tmap.NumMyElements();
00462 
00463     // importer
00464     Epetra_Import importer(tmap, smap);
00465 
00466     // number of non zeros in my new rows
00467     Epetra_Vector numOld(smap);
00468     Epetra_Vector numNew(tmap);
00469 
00470     for (int lid=0; lid < oldNumElts; lid++){
00471       numOld[lid] = from->NumMyIndices(lid);
00472     }
00473     numNew.Import(numOld, importer, Insert);
00474 
00475     Array<int> nnz(newNumElts);
00476     for (int lid=0; lid < newNumElts; lid++){
00477       nnz[lid] = static_cast<int>(numNew[lid]);
00478     }
00479 
00480     // target graph
00481     RCP<Epetra_CrsGraph> G = rcp(
00482       new Epetra_CrsGraph(Copy, tmap, nnz.getRawPtr(), true));
00483     G->Import(*from, importer, Insert);
00484     G->FillComplete();
00485 
00486     return G;
00487   }
00488 
00489 };
00490 
00492 // Xpetra::CrsGraph
00493 // KDDKDD Do we need specializations for Xpetra::TpetraCrsGraph and
00494 // KDDKDD Xpetra::EpetraCrsGraph?
00495 template <typename lno_t,
00496           typename gno_t,
00497           typename node_t>
00498 struct XpetraTraits<Xpetra::CrsGraph<lno_t, gno_t, node_t> >
00499 {
00500   typedef Xpetra::CrsGraph<lno_t, gno_t, node_t> x_graph_t;
00501   typedef Xpetra::TpetraCrsGraph<lno_t, gno_t, node_t> xt_graph_t;
00502   typedef Tpetra::CrsGraph<lno_t,gno_t,node_t> t_graph_t; 
00503 
00504   static inline RCP<const x_graph_t>
00505     convertToXpetra(const RCP<const x_graph_t> &a)
00506     {
00507       return a;
00508     }
00509 
00510   static RCP<const x_graph_t> doMigration(const RCP<const x_graph_t> &from,
00511       lno_t numLocalRows, const gno_t *myNewRows)
00512   {
00513     Xpetra::UnderlyingLib lib = from->getRowMap()->lib();
00514 
00515     if (lib == Xpetra::UseEpetra){
00516        throw logic_error("compiler should have used specialization");
00517     } else{
00518       // Do the import with the Tpetra::CrsGraph traits object
00519       const x_graph_t *xg = from.get();
00520       const xt_graph_t *xtg = dynamic_cast<const xt_graph_t *>(xg);
00521       RCP<const t_graph_t> tg = xtg->getTpetra_CrsGraph();
00522 
00523       RCP<const t_graph_t> tgnew = XpetraTraits<t_graph_t>::doMigration(
00524         tg, numLocalRows, myNewRows);
00525 
00526       RCP<const x_graph_t> xgnew =
00527         XpetraTraits<t_graph_t>::convertToXpetra(tgnew);
00528       return xgnew;
00529     }
00530   }
00531 };
00532 
00533 
00534 
00536 // Xpetra::RowMatrix
00537 template <typename scalar_t,
00538           typename lno_t,
00539           typename gno_t,
00540           typename node_t>
00541 struct XpetraTraits<Xpetra::RowMatrix<scalar_t, lno_t, gno_t, node_t> >
00542 {
00543   typedef Xpetra::RowMatrix<scalar_t, lno_t, gno_t, node_t> x_matrix_t;
00544   typedef Xpetra::TpetraRowMatrix<scalar_t, lno_t, gno_t, node_t> xt_matrix_t;
00545   typedef Tpetra::RowMatrix<scalar_t,lno_t,gno_t,node_t> t_matrix_t; 
00546 
00547   static inline RCP<const x_matrix_t>
00548     convertToXpetra( const RCP<const x_matrix_t > &a)
00549     {
00550       return a;
00551     }
00552 
00553   static RCP<const x_matrix_t> doMigration(const RCP<const x_matrix_t> &from,
00554       lno_t numLocalRows, const gno_t *myNewRows)
00555   {
00556     Xpetra::UnderlyingLib lib = from->getRowMap()->lib();
00557 
00558     if (lib == Xpetra::UseEpetra){
00559        throw logic_error("compiler should have used specialization");
00560     } else{
00561       // Do the import with the Tpetra::CrsMatrix traits object
00562       const x_matrix_t *xm = from.get();
00563       const xt_matrix_t *xtm = dynamic_cast<const xt_matrix_t *>(xm);
00564       RCP<const t_matrix_t> tm = xtm->getTpetra_CrsMatrix();
00565 
00566       RCP<const t_matrix_t> tmnew = XpetraTraits<t_matrix_t>::doMigration(
00567         tm, numLocalRows, myNewRows);
00568 
00569       RCP<const x_matrix_t> xmnew = 
00570         XpetraTraits<t_matrix_t>::convertToXpetra(tmnew);
00571 
00572       return xmnew;
00573     }
00574   }
00575 };
00576 
00578 // Xpetra::CrsGraph specialization
00579 template < typename node_t>
00580 struct XpetraTraits<Xpetra::CrsGraph<int, int, node_t> >
00581 {
00582   typedef int lno_t;
00583   typedef int gno_t;
00584   typedef Xpetra::CrsGraph<lno_t, gno_t, node_t> x_graph_t;
00585   typedef Xpetra::TpetraCrsGraph<lno_t, gno_t, node_t> xt_graph_t;
00586   typedef Tpetra::CrsGraph<lno_t,gno_t,node_t> t_graph_t; 
00587   typedef Xpetra::EpetraCrsGraph xe_graph_t;
00588   typedef Epetra_CrsGraph e_graph_t; 
00589 
00590   static inline RCP<const x_graph_t>
00591     convertToXpetra(const RCP<const x_graph_t> &a)
00592     {
00593       return a;
00594     }
00595 
00596   static RCP<const x_graph_t> doMigration(const RCP<const x_graph_t> &from,
00597       lno_t numLocalRows, const gno_t *myNewRows)
00598   {
00599     Xpetra::UnderlyingLib lib = from->getRowMap()->lib();
00600     const x_graph_t *xg = from.get();
00601 
00602     if (lib == Xpetra::UseEpetra){
00603       // Do the import with the Epetra_CrsGraph traits object
00604       const xe_graph_t *xeg = dynamic_cast<const xe_graph_t *>(xg);
00605       RCP<const e_graph_t> eg = xeg->getEpetra_CrsGraph();
00606 
00607       RCP<const e_graph_t> egnew = XpetraTraits<e_graph_t>::doMigration(
00608         eg, numLocalRows, myNewRows);
00609 
00610       RCP<const x_graph_t> xgnew =
00611         XpetraTraits<e_graph_t>::convertToXpetra(egnew);
00612 
00613       return xgnew;
00614     } else{
00615       // Do the import with the Tpetra::CrsGraph traits object
00616       const xt_graph_t *xtg = dynamic_cast<const xt_graph_t *>(xg);
00617       RCP<const t_graph_t> tg = xtg->getTpetra_CrsGraph();
00618 
00619       RCP<const t_graph_t> tgnew = XpetraTraits<t_graph_t>::doMigration(
00620         tg, numLocalRows, myNewRows);
00621 
00622       RCP<const x_graph_t> xgnew =
00623         XpetraTraits<t_graph_t>::convertToXpetra(tgnew);
00624 
00625       return xgnew;
00626     }
00627   }
00628 };
00629 
00631 // Tpetra::Vector
00632 template <typename scalar_t,
00633           typename lno_t,
00634           typename gno_t,
00635           typename node_t>
00636 struct XpetraTraits<Tpetra::Vector<scalar_t, lno_t, gno_t, node_t> >
00637 {
00638   typedef Tpetra::Vector<scalar_t, lno_t, gno_t, node_t> t_vector_t;
00639   typedef Xpetra::TpetraVector<scalar_t, lno_t, gno_t, node_t> xt_vector_t;
00640   typedef Xpetra::Vector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
00641 
00642   static inline RCP<const x_vector_t>
00643     convertToXpetra(const RCP<const t_vector_t> &a)
00644     {
00645       return rcp(new xt_vector_t(rcp_const_cast<t_vector_t>(a)));
00646     }
00647 
00648   static RCP<const t_vector_t> doMigration(const RCP<const t_vector_t> &from,
00649       lno_t numLocalElts, const gno_t *myNewElts)
00650   {
00651     lno_t base = 0;
00652     typedef Tpetra::Map<lno_t, gno_t, node_t> map_t;
00653 
00654     // source map
00655     const RCP<const map_t> &smap = from->getMap();
00656     gno_t numGlobalElts = smap->getGlobalNumElements();
00657 
00658     // target map
00659     ArrayView<const gno_t> eltList(myNewElts, numLocalElts);
00660     const RCP<const Teuchos::Comm<int> > comm = from->getMap()->getComm();
00661     RCP<const map_t> tmap = rcp(
00662       new map_t(numGlobalElts, eltList, base, comm));
00663 
00664     // importer
00665     Tpetra::Import<lno_t, gno_t, node_t> importer(smap, tmap);
00666 
00667     // target vector 
00668     RCP<t_vector_t> V = 
00669       Tpetra::createVector<scalar_t,lno_t,gno_t,node_t>(tmap);
00670     V->doImport(*from, importer, Tpetra::INSERT);
00671 
00672     return V;
00673   }
00674 };
00675 
00677 // Epetra_Vector
00678 template < >
00679 struct XpetraTraits<Epetra_Vector>
00680 {
00681   typedef InputTraits<Epetra_Vector>::lno_t    lno_t;
00682   typedef InputTraits<Epetra_Vector>::gno_t    gno_t;
00683   typedef InputTraits<Epetra_Vector>::node_t   node_t;
00684   typedef InputTraits<Epetra_Vector>::scalar_t   scalar_t;
00685   
00686   typedef Xpetra::Vector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
00687 
00688   static inline RCP<const x_vector_t>
00689     convertToXpetra(const RCP<const Epetra_Vector> &a)
00690     {
00691       RCP<Xpetra::EpetraVector> xev = rcp(
00692         new Xpetra::EpetraVector(rcp_const_cast<Epetra_Vector>(a)));
00693       return rcp_implicit_cast<x_vector_t>(xev);
00694     }
00695 
00696   static RCP<Epetra_Vector> doMigration(const RCP<const Epetra_Vector> &from,
00697       lno_t numLocalElts, const gno_t *myNewElts)
00698   {
00699     lno_t base = 0;
00700     // source map
00701     const Epetra_BlockMap &smap = from->Map();
00702     gno_t numGlobalElts = smap.NumGlobalElements();
00703 
00704     // target map
00705     const Epetra_Comm &comm = from->Comm();
00706     const Epetra_BlockMap tmap(numGlobalElts, numLocalElts, myNewElts, 
00707       1, base, comm);
00708 
00709     // importer
00710     Epetra_Import importer(tmap, smap);
00711 
00712     // target vector 
00713     RCP<Epetra_Vector> V = rcp(new Epetra_Vector(tmap, true));
00714     Epetra_CombineMode c = Insert;
00715     V->Import(*from, importer, c);
00716 
00717     return V;
00718   }
00719 };
00720 
00722 // Xpetra::Vector
00723 template <typename scalar_t,
00724           typename lno_t,
00725           typename gno_t,
00726           typename node_t>
00727 struct XpetraTraits<Xpetra::Vector<scalar_t, lno_t, gno_t, node_t> >
00728 {
00729   typedef Xpetra::Vector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
00730   typedef Xpetra::TpetraVector<scalar_t, lno_t, gno_t, node_t> xt_vector_t;
00731   typedef Tpetra::Vector<scalar_t, lno_t, gno_t, node_t> t_vector_t;
00732 
00733   static inline RCP<const x_vector_t>
00734     convertToXpetra(const RCP<const x_vector_t> &a)
00735     {
00736       return a;
00737     }
00738 
00739   static RCP<const x_vector_t> doMigration(const RCP<const x_vector_t> &from,
00740       lno_t numLocalRows, const gno_t *myNewRows)
00741   {
00742     Xpetra::UnderlyingLib lib = from->getMap()->lib();
00743 
00744     if (lib == Xpetra::UseEpetra){
00745        throw logic_error("compiler should have used specialization");
00746     } else{
00747       // Do the import with the Tpetra::Vector traits object
00748       const x_vector_t *xv = from.get();
00749       const xt_vector_t *xtv = dynamic_cast<const xt_vector_t *>(xv);
00750       RCP<const t_vector_t> tv = xtv->getTpetra_Vector();
00751 
00752       RCP<const t_vector_t> tvnew = XpetraTraits<t_vector_t>::doMigration(
00753         tv, numLocalRows, myNewRows);
00754 
00755       RCP<const x_vector_t> xvnew =
00756         XpetraTraits<t_vector_t>::convertToXpetra(tvnew);
00757 
00758       return xvnew;
00759     }
00760   }
00761 };
00762 
00764 // Xpetra::Vector specialization
00765 template <typename node_t>
00766 struct XpetraTraits<Xpetra::Vector<double, int, int, node_t> >
00767 {
00768   typedef double scalar_t;
00769   typedef int lno_t;
00770   typedef int gno_t;
00771   typedef Xpetra::Vector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
00772   typedef Xpetra::TpetraVector<scalar_t, lno_t, gno_t, node_t> xt_vector_t;
00773   typedef Tpetra::Vector<scalar_t, lno_t, gno_t, node_t> t_vector_t;
00774   typedef Xpetra::EpetraVector xe_vector_t;
00775   typedef Epetra_Vector e_vector_t;
00776 
00777   static inline RCP<const x_vector_t>
00778     convertToXpetra(const RCP<const x_vector_t> &a)
00779     {
00780       return a;
00781     }
00782 
00783   static RCP<const x_vector_t> doMigration(const RCP<const x_vector_t> &from,
00784       lno_t numLocalRows, const gno_t *myNewRows)
00785   {
00786     Xpetra::UnderlyingLib lib = from->getMap()->lib();
00787     const x_vector_t *vec = from.get();
00788 
00789     if (lib == Xpetra::UseEpetra){
00790       // Do the import with the Epetra_Vector traits object
00791       const xe_vector_t *xev = dynamic_cast<const xe_vector_t *>(vec);
00792       RCP<const e_vector_t> ev = rcp(xev->getEpetra_Vector());
00793 
00794       RCP<const e_vector_t> evnew = XpetraTraits<e_vector_t>::doMigration(
00795         ev, numLocalRows, myNewRows);
00796 
00797       RCP<const x_vector_t> xvnew =
00798         XpetraTraits<e_vector_t>::convertToXpetra(evnew);
00799           
00800       return xvnew;
00801     } else{
00802       // Do the import with the Tpetra::Vector traits object
00803       const xt_vector_t *xtv = dynamic_cast<const xt_vector_t *>(vec);
00804       RCP<t_vector_t> tv = xtv->getTpetra_Vector();
00805       RCP<const t_vector_t> ctv = rcp_const_cast<const t_vector_t>(tv);
00806 
00807       RCP<const t_vector_t> tvnew = XpetraTraits<t_vector_t>::doMigration(
00808         ctv, numLocalRows, myNewRows);
00809 
00810       RCP<const x_vector_t> xvnew =
00811         XpetraTraits<t_vector_t>::convertToXpetra(tvnew);
00812 
00813       return xvnew;
00814     }
00815   }
00816 };
00817 
00819 // Tpetra::MultiVector
00820 template <typename scalar_t,
00821           typename lno_t,
00822           typename gno_t,
00823           typename node_t>
00824 struct XpetraTraits<Tpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> >
00825 {
00826   typedef Tpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> t_vector_t;
00827   typedef Xpetra::TpetraMultiVector<scalar_t, lno_t, gno_t, node_t> xt_vector_t;
00828   typedef Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
00829 
00830   static inline RCP<const x_vector_t>
00831     convertToXpetra(const RCP<const t_vector_t> &a)
00832     {
00833       return rcp(new xt_vector_t(rcp_const_cast<t_vector_t>(a)));
00834     }
00835 
00836   static RCP<const t_vector_t> doMigration(const RCP<const t_vector_t> &from,
00837       lno_t numLocalElts, const gno_t *myNewElts)
00838   {
00839     typedef Tpetra::Map<lno_t, gno_t, node_t> map_t;
00840     lno_t base = 0;
00841 
00842     // source map
00843     const RCP<const map_t> &smap = from->getMap();
00844     gno_t numGlobalElts = smap->getGlobalNumElements();
00845 
00846     // target map
00847     ArrayView<const gno_t> eltList(myNewElts, numLocalElts);
00848     const RCP<const Teuchos::Comm<int> > comm = from->getMap()->getComm();
00849     RCP<const map_t> tmap = rcp(
00850       new map_t(numGlobalElts, eltList, base, comm));
00851 
00852     // importer
00853     Tpetra::Import<lno_t, gno_t, node_t> importer(smap, tmap);
00854 
00855     // target vector 
00856     RCP<t_vector_t> MV = rcp(
00857       new t_vector_t(tmap, from->getNumVectors(), true));
00858     MV->doImport(*from, importer, Tpetra::INSERT);
00859 
00860     return MV;
00861   }
00862 };
00863 
00865 // Epetra_MultiVector
00866 template < >
00867 struct XpetraTraits<Epetra_MultiVector>
00868 {
00869   typedef InputTraits<Epetra_MultiVector>::lno_t    lno_t;
00870   typedef InputTraits<Epetra_MultiVector>::gno_t    gno_t;
00871   typedef InputTraits<Epetra_MultiVector>::node_t   node_t;
00872   typedef InputTraits<Epetra_MultiVector>::scalar_t   scalar_t;
00873   typedef Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> x_mvector_t;
00874 
00875   static inline RCP<const x_mvector_t>
00876     convertToXpetra(const RCP<const Epetra_MultiVector> &a)
00877     {
00878       RCP<Xpetra::EpetraMultiVector> xemv = rcp(
00879         new Xpetra::EpetraMultiVector(
00880           rcp_const_cast<Epetra_MultiVector>(a)));
00881       return rcp_implicit_cast<x_mvector_t>(xemv);
00882     }
00883 
00884   static RCP<Epetra_MultiVector> doMigration(
00885     const RCP<const Epetra_MultiVector> &from,
00886     lno_t numLocalElts, const gno_t *myNewElts)
00887   {
00888     lno_t base = 0;
00889     // source map
00890     const Epetra_BlockMap &smap = from->Map();
00891     gno_t numGlobalElts = smap.NumGlobalElements();
00892 
00893     // target map
00894     const Epetra_Comm &comm = from->Comm();
00895     const Epetra_BlockMap tmap(numGlobalElts, numLocalElts, myNewElts, 
00896       1, base, comm);
00897 
00898     // importer
00899     Epetra_Import importer(tmap, smap);
00900 
00901     // target vector 
00902     RCP<Epetra_MultiVector> MV = rcp(
00903       new Epetra_MultiVector(tmap, from->NumVectors(), true));
00904     Epetra_CombineMode c = Insert;
00905     MV->Import(*from, importer, c);
00906 
00907     return MV;
00908   }
00909 };
00910 
00912 // Xpetra::MultiVector
00913 template <typename scalar_t,
00914           typename lno_t,
00915           typename gno_t,
00916           typename node_t>
00917 struct XpetraTraits<Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> >
00918 {
00919   typedef Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> x_mvector_t;
00920   typedef Xpetra::TpetraMultiVector<scalar_t, lno_t, gno_t, node_t> xt_mvector_t;
00921   typedef Tpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> t_mvector_t;
00922 
00923   static inline RCP<const x_mvector_t>
00924     convertToXpetra(const RCP<const x_mvector_t> &a)
00925     {
00926       return a;
00927     }
00928 
00929   static RCP<const x_mvector_t> doMigration(const RCP<const x_mvector_t> &from,
00930       lno_t numLocalRows, const gno_t *myNewRows)
00931   {
00932     Xpetra::UnderlyingLib lib = from->getMap()->lib();
00933 
00934     if (lib == Xpetra::UseEpetra){
00935        throw logic_error("compiler should have used specialization");
00936     } else{
00937       // Do the import with the Tpetra::MultiVector traits object
00938       const x_mvector_t *xmv = from.get();
00939       const xt_mvector_t *xtv = dynamic_cast<const xt_mvector_t *>(xmv);
00940       RCP<t_mvector_t> tv = xtv->getTpetra_MultiVector();
00941       RCP<const t_mvector_t> ctv = rcp_const_cast<const t_mvector_t>(tv);
00942 
00943       RCP<const t_mvector_t> tvnew = XpetraTraits<t_mvector_t>::doMigration(
00944         ctv, numLocalRows, myNewRows);
00945 
00946       RCP<const x_mvector_t> xvnew =
00947         XpetraTraits<t_mvector_t>::convertToXpetra(tvnew);
00948 
00949       return xvnew;
00950     }
00951   }
00952 };
00953 
00955 // Xpetra::MultiVector specialization
00956 template <typename node_t>
00957 struct XpetraTraits<Xpetra::MultiVector<double, int, int, node_t> >
00958 {
00959   typedef double scalar_t;
00960   typedef int lno_t;
00961   typedef int gno_t;
00962   typedef Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> x_mvector_t;
00963   typedef Xpetra::TpetraMultiVector<scalar_t, lno_t, gno_t, node_t> xt_mvector_t;
00964   typedef Tpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> t_mvector_t;
00965   typedef Xpetra::EpetraMultiVector xe_mvector_t;
00966   typedef Epetra_MultiVector e_mvector_t;
00967 
00968   static inline RCP<const x_mvector_t>
00969     convertToXpetra(const RCP<const x_mvector_t> &a)
00970     {
00971       return a;
00972     }
00973 
00974   static RCP<const x_mvector_t> doMigration(const RCP<const x_mvector_t> &from,
00975       lno_t numLocalRows, const gno_t *myNewRows)
00976   {
00977     Xpetra::UnderlyingLib lib = from->getMap()->lib();
00978     const x_mvector_t *xmv = from.get();
00979 
00980     if (lib == Xpetra::UseEpetra){
00981       // Do the import with the Epetra_MultiVector traits object
00982       const xe_mvector_t *xev = dynamic_cast<const xe_mvector_t *>(xmv);
00983       RCP<e_mvector_t> ev = xev->getEpetra_MultiVector();
00984       RCP<const e_mvector_t> cev = rcp_const_cast<const e_mvector_t>(ev);
00985 
00986       RCP<const e_mvector_t> evnew = XpetraTraits<e_mvector_t>::doMigration(
00987         cev, numLocalRows, myNewRows);
00988 
00989       RCP<const x_mvector_t> xvnew =
00990         XpetraTraits<e_mvector_t>::convertToXpetra(evnew);
00991 
00992       return xvnew;
00993 
00994     } else{
00995       // Do the import with the Tpetra::MultiVector traits object
00996       const xt_mvector_t *xtv = dynamic_cast<const xt_mvector_t *>(xmv);
00997       RCP<t_mvector_t> tv = xtv->getTpetra_MultiVector();
00998       RCP<const t_mvector_t> ctv = rcp_const_cast<const t_mvector_t>(tv);
00999 
01000       RCP<const t_mvector_t> tvnew = XpetraTraits<t_mvector_t>::doMigration(
01001         ctv, numLocalRows, myNewRows);
01002 
01003       RCP<const x_mvector_t> xvnew =
01004         XpetraTraits<t_mvector_t>::convertToXpetra(tvnew);
01005 
01006       return xvnew;
01007     }
01008   }
01009 };
01010 
01011 #endif // DOXYGEN_SHOULD_SKIP_THIS
01012 
01013 }  //namespace Zoltan2
01014 
01015 #endif // _ZOLTAN2_XPETRATRAITS_HPP_