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