50 #ifndef _ZOLTAN2_XPETRATRAITS_HPP_
51 #define _ZOLTAN2_XPETRATRAITS_HPP_
56 #include <Xpetra_TpetraCrsGraph.hpp>
57 #include <Xpetra_TpetraCrsMatrix.hpp>
58 #include <Xpetra_TpetraVector.hpp>
59 #include <Tpetra_Vector.hpp>
61 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA)
62 #include <Xpetra_EpetraCrsMatrix.hpp>
63 #include <Xpetra_EpetraVector.hpp>
64 #include <Xpetra_EpetraUtils.hpp>
94 template <
typename User>
118 size_t numLocalRows,
const gno_t *myNewRows)
120 return Teuchos::null;
124 #ifndef DOXYGEN_SHOULD_SKIP_THIS
128 template <
typename scalar_t,
132 struct XpetraTraits<Tpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> >
134 typedef typename Xpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t>
xmatrix_t;
135 typedef typename Xpetra::TpetraCrsMatrix<scalar_t,lno_t,gno_t,node_t> xtmatrix_t;
136 typedef typename Tpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t>
tmatrix_t;
140 return rcp(
new xtmatrix_t(a));
144 size_t numLocalRows,
const gno_t *myNewRows)
146 typedef Tpetra::Map<lno_t, gno_t, node_t>
map_t;
149 const RCP<const map_t> &smap = from.getRowMap();
150 gno_t numGlobalRows = smap->getGlobalNumElements();
151 gno_t base = smap->getMinAllGlobalIndex();
154 ArrayView<const gno_t> rowList(myNewRows, numLocalRows);
155 const RCP<const Teuchos::Comm<int> > &comm = from.getComm();
156 RCP<const map_t> tmap = rcp(
new map_t(numGlobalRows, rowList, base, comm));
159 Tpetra::Import<lno_t, gno_t, node_t> importer(smap, tmap);
176 from.importAndFillComplete(M, importer, tmap, tmap);
212 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA)
215 struct XpetraTraits<Epetra_CrsMatrix>
217 typedef InputTraits<Epetra_CrsMatrix>::scalar_t scalar_t;
222 static inline RCP<Xpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t> >
225 RCP<Xpetra::EpetraCrsMatrixT<gno_t, node_t> > xa;
227 xa = rcp(
new Xpetra::EpetraCrsMatrixT<gno_t, node_t>(a));
229 catch (std::exception &e) {
230 if (std::is_same<node_t, Xpetra::EpetraNode>::value)
231 throw std::runtime_error(std::string(
"Cannot convert from "
232 "Epetra_CrsMatrix to "
233 "Xpetra::EpetraCrsMatrixT\n")
236 throw std::runtime_error(std::string(
"Cannot convert from "
237 "Epetra_CrsMatrix to "
238 "Xpetra::EpetraCrsMatrixT\n"
239 "Use node_t that is supported by "
240 "Xpetra with Epetra classes\n")
247 static RCP<Epetra_CrsMatrix>
doMigration(
const Epetra_CrsMatrix &from,
248 size_t numLocalRows,
const gno_t *myNewRows)
251 const Epetra_Map &smap = from.RowMap();
252 gno_t numGlobalRows = smap.NumGlobalElements();
253 int base = smap.MinAllGID();
256 const Epetra_Comm &comm = from.Comm();
257 Epetra_Map tmap(numGlobalRows, numLocalRows, myNewRows, base, comm);
260 Epetra_Import importer(tmap, smap);
272 RCP<Epetra_CrsMatrix> M = rcp(
new Epetra_CrsMatrix(from, importer,
308 template <
typename scalar_t,
312 struct XpetraTraits<Xpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> >
314 typedef Xpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> x_matrix_t;
315 typedef Xpetra::TpetraCrsMatrix<scalar_t, lno_t, gno_t, node_t> xt_matrix_t;
316 typedef Tpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t> t_matrix_t;
318 static inline RCP<x_matrix_t>
convertToXpetra(
const RCP<x_matrix_t > &a)
323 static RCP<x_matrix_t>
doMigration(
const x_matrix_t &from,
324 size_t numLocalRows,
const gno_t *myNewRows)
326 Xpetra::UnderlyingLib lib = from.getRowMap()->lib();
328 if (lib == Xpetra::UseEpetra){
329 throw std::logic_error(
"compiler should have used specialization");
332 const xt_matrix_t *xtm =
dynamic_cast<const xt_matrix_t *
>(&from);
333 RCP<const t_matrix_t> tm = xtm->getTpetra_CrsMatrix();
336 *tm, numLocalRows, myNewRows);
348 template <
typename node_t>
349 struct XpetraTraits<Xpetra::CrsMatrix<double, int, int, node_t> >
351 typedef double scalar_t;
354 typedef Xpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> x_matrix_t;
355 typedef Xpetra::TpetraCrsMatrix<scalar_t, lno_t, gno_t, node_t> xt_matrix_t;
356 typedef Tpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t> t_matrix_t;
358 static inline RCP<x_matrix_t>
convertToXpetra(
const RCP<x_matrix_t > &a)
363 static RCP<x_matrix_t>
doMigration(
const x_matrix_t &from,
364 size_t numLocalRows,
const gno_t *myNewRows)
366 Xpetra::UnderlyingLib lib = from.getRowMap()->lib();
368 if (lib == Xpetra::UseEpetra){
369 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA)
370 typedef Epetra_CrsMatrix e_matrix_t;
371 typedef Xpetra::EpetraCrsMatrixT<gno_t,node_t> xe_matrix_t;
373 const xe_matrix_t *xem =
dynamic_cast<const xe_matrix_t *
>(&from);
374 RCP<const e_matrix_t> em = xem->getEpetra_CrsMatrix();
377 *em, numLocalRows, myNewRows);
383 throw std::runtime_error(
"Xpetra with Epetra requested, but "
384 "Trilinos is not built with Epetra");
388 const xt_matrix_t *xtm =
dynamic_cast<const xt_matrix_t *
>(&from);
389 RCP<const t_matrix_t> tm = xtm->getTpetra_CrsMatrix();
392 *tm, numLocalRows, myNewRows);
404 template <
typename lno_t,
407 struct XpetraTraits<Tpetra::CrsGraph<lno_t, gno_t, node_t> >
409 typedef typename Xpetra::CrsGraph<lno_t, gno_t, node_t>
xgraph_t;
410 typedef typename Xpetra::TpetraCrsGraph<lno_t, gno_t, node_t> xtgraph_t;
411 typedef typename Tpetra::CrsGraph<lno_t, gno_t, node_t>
tgraph_t;
415 return rcp(
new xtgraph_t(a));
419 size_t numLocalRows,
const gno_t *myNewRows)
421 typedef Tpetra::Map<lno_t, gno_t, node_t>
map_t;
424 const RCP<const map_t> &smap = from.getRowMap();
425 int oldNumElts = smap->getLocalNumElements();
426 gno_t numGlobalRows = smap->getGlobalNumElements();
427 gno_t base = smap->getMinAllGlobalIndex();
430 ArrayView<const gno_t> rowList(myNewRows, numLocalRows);
431 const RCP<const Teuchos::Comm<int> > &comm = from.getComm();
432 RCP<const map_t> tmap = rcp(
new map_t(numGlobalRows, rowList, base, comm));
435 Tpetra::Import<lno_t, gno_t, node_t> importer(smap, tmap);
438 typedef Tpetra::Vector<gno_t, lno_t, gno_t, node_t> vector_t;
439 vector_t numOld(smap);
440 vector_t numNew(tmap);
441 for (
int lid=0; lid < oldNumElts; lid++){
442 numOld.replaceGlobalValue(smap->getGlobalElement(lid),
443 from.getNumEntriesInLocalRow(lid));
445 numNew.doImport(numOld, importer, Tpetra::INSERT);
447 size_t numElts = tmap->getLocalNumElements();
448 ArrayRCP<const gno_t> nnz;
450 nnz = numNew.getData(0);
452 ArrayRCP<const size_t> nnz_size_t;
454 if (numElts &&
sizeof(
gno_t) !=
sizeof(
size_t)){
455 size_t *vals =
new size_t [numElts];
456 nnz_size_t = arcp(vals, 0, numElts,
true);
457 for (
size_t i=0; i < numElts; i++){
458 vals[i] =
static_cast<size_t>(nnz[i]);
462 nnz_size_t = arcp_reinterpret_cast<
const size_t>(nnz);
466 RCP<tgraph_t> G = rcp(
new tgraph_t(tmap, nnz_size_t()));
468 G->doImport(from, importer, Tpetra::INSERT);
477 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA)
480 struct XpetraTraits<Epetra_CrsGraph>
485 static inline RCP<Xpetra::CrsGraph<lno_t,gno_t,node_t> >
488 RCP<Xpetra::EpetraCrsGraphT<gno_t, node_t> > xa;
490 xa = rcp(
new Xpetra::EpetraCrsGraphT<gno_t, node_t>(a));
492 catch (std::exception &e) {
493 if (std::is_same<node_t, Xpetra::EpetraNode>::value)
494 throw std::runtime_error(std::string(
"Cannot convert from "
495 "Epetra_CrsGraph to "
496 "Xpetra::EpetraCrsGraphT\n")
499 throw std::runtime_error(std::string(
"Cannot convert from "
500 "Epetra_CrsGraph to "
501 "Xpetra::EpetraCrsGraphT\n"
502 "Use node_t that is supported by "
503 "Xpetra with Epetra classes\n")
509 static RCP<Epetra_CrsGraph>
doMigration(
const Epetra_CrsGraph &from,
510 size_t numLocalRows,
const gno_t *myNewRows)
513 const Epetra_BlockMap &smap = from.RowMap();
514 gno_t numGlobalRows = smap.NumGlobalElements();
515 lno_t oldNumElts = smap.NumMyElements();
516 int base = smap.MinAllGID();
519 const Epetra_Comm &comm = from.Comm();
520 Epetra_BlockMap tmap(numGlobalRows, numLocalRows, myNewRows, 1, base, comm);
521 lno_t newNumElts = tmap.NumMyElements();
524 Epetra_Import importer(tmap, smap);
527 Epetra_Vector numOld(smap);
528 Epetra_Vector numNew(tmap);
530 for (
int lid=0; lid < oldNumElts; lid++){
531 numOld[lid] = from.NumMyIndices(lid);
533 numNew.Import(numOld, importer, Insert);
535 Array<int> nnz(newNumElts);
536 for (
int lid=0; lid < newNumElts; lid++){
537 nnz[lid] =
static_cast<int>(numNew[lid]);
541 RCP<Epetra_CrsGraph> G = rcp(
new Epetra_CrsGraph(::Copy, tmap, nnz.getRawPtr(),
true));
542 G->Import(from, importer, Insert);
555 template <
typename lno_t,
558 struct XpetraTraits<Xpetra::CrsGraph<lno_t, gno_t, node_t> >
560 typedef Xpetra::CrsGraph<lno_t, gno_t, node_t> x_graph_t;
561 typedef Xpetra::TpetraCrsGraph<lno_t, gno_t, node_t> xt_graph_t;
562 typedef Tpetra::CrsGraph<lno_t,gno_t,node_t> t_graph_t;
569 static RCP<x_graph_t>
doMigration(
const x_graph_t &from,
570 size_t numLocalRows,
const gno_t *myNewRows)
572 Xpetra::UnderlyingLib lib = from.getRowMap()->lib();
574 if (lib == Xpetra::UseEpetra){
575 throw std::logic_error(
"compiler should have used specialization");
578 const xt_graph_t *xtg =
dynamic_cast<const xt_graph_t *
>(&from);
579 RCP<const t_graph_t> tg = xtg->getTpetra_CrsGraph();
582 *tg, numLocalRows, myNewRows);
593 template <
typename node_t>
594 struct XpetraTraits<Xpetra::CrsGraph<int, int, node_t> >
598 typedef Xpetra::CrsGraph<lno_t, gno_t, node_t> x_graph_t;
599 typedef Xpetra::TpetraCrsGraph<lno_t, gno_t, node_t> xt_graph_t;
600 typedef Tpetra::CrsGraph<lno_t,gno_t,node_t> t_graph_t;
607 static RCP<x_graph_t>
doMigration(
const x_graph_t &from,
608 size_t numLocalRows,
const gno_t *myNewRows)
610 Xpetra::UnderlyingLib lib = from.getRowMap()->lib();
612 if (lib == Xpetra::UseEpetra){
613 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA)
614 typedef Xpetra::EpetraCrsGraphT<gno_t,node_t> xe_graph_t;
615 typedef Epetra_CrsGraph e_graph_t;
617 const xe_graph_t *xeg =
dynamic_cast<const xe_graph_t *
>(&from);
618 RCP<const e_graph_t> eg = xeg->getEpetra_CrsGraph();
621 *eg, numLocalRows, myNewRows);
627 throw std::runtime_error(
"Xpetra with Epetra requested, but "
628 "Trilinos is not built with Epetra");
632 const xt_graph_t *xtg =
dynamic_cast<const xt_graph_t *
>(&from);
633 RCP<const t_graph_t> tg = xtg->getTpetra_CrsGraph();
636 *tg, numLocalRows, myNewRows);
647 template <
typename scalar_t,
651 struct XpetraTraits<Tpetra::
Vector<scalar_t, lno_t, gno_t, node_t> >
653 typedef Tpetra::Vector<scalar_t, lno_t, gno_t, node_t> t_vector_t;
654 typedef Xpetra::TpetraVector<scalar_t, lno_t, gno_t, node_t> xt_vector_t;
655 typedef Xpetra::Vector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
657 static inline RCP<x_vector_t>
convertToXpetra(
const RCP<t_vector_t> &a)
659 return rcp(
new xt_vector_t(a));
662 static RCP<t_vector_t>
doMigration(
const t_vector_t &from,
663 size_t numLocalElts,
const gno_t *myNewElts)
665 typedef Tpetra::Map<lno_t, gno_t, node_t>
map_t;
668 const RCP<const map_t> &smap = from.getMap();
669 gno_t numGlobalElts = smap->getGlobalNumElements();
670 gno_t base = smap->getMinAllGlobalIndex();
673 ArrayView<const gno_t> eltList(myNewElts, numLocalElts);
674 const RCP<const Teuchos::Comm<int> > comm = from.getMap()->getComm();
675 RCP<const map_t> tmap = rcp(
new map_t(numGlobalElts, eltList, base, comm));
678 Tpetra::Import<lno_t, gno_t, node_t> importer(smap, tmap);
682 Tpetra::createVector<scalar_t,lno_t,gno_t,node_t>(tmap);
683 V->doImport(from, importer, Tpetra::INSERT);
690 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA)
693 struct XpetraTraits<Epetra_Vector>
698 typedef InputTraits<Epetra_Vector>::scalar_t scalar_t;
700 typedef Xpetra::Vector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
702 static inline RCP<x_vector_t>
convertToXpetra(
const RCP<Epetra_Vector> &a)
704 RCP<Xpetra::EpetraVectorT<gno_t, node_t> > xev;
706 xev = rcp(
new Xpetra::EpetraVectorT<gno_t,node_t>(a));
708 catch (std::exception &e) {
709 if (std::is_same<node_t, Xpetra::EpetraNode>::value)
710 throw std::runtime_error(std::string(
"Cannot convert from "
712 "Xpetra::EpetraVectorT\n")
715 throw std::runtime_error(std::string(
"Cannot convert from "
717 "Xpetra::EpetraVectorT\n"
718 "Use node_t that is supported by "
719 "Xpetra with Epetra classes\n")
722 return rcp_implicit_cast<x_vector_t>(xev);
725 static RCP<Epetra_Vector>
doMigration(
const Epetra_Vector &from,
726 size_t numLocalElts,
const gno_t *myNewElts)
729 const Epetra_BlockMap &smap = from.Map();
730 gno_t numGlobalElts = smap.NumGlobalElements();
731 int base = smap.MinAllGID();
734 const Epetra_Comm &comm = from.Comm();
735 const Epetra_BlockMap tmap(numGlobalElts, numLocalElts, myNewElts,
739 Epetra_Import importer(tmap, smap);
742 RCP<Epetra_Vector> V = rcp(
new Epetra_Vector(tmap,
true));
743 Epetra_CombineMode c = Insert;
744 V->Import(from, importer, c);
753 template <
typename scalar_t,
757 struct XpetraTraits<Xpetra::
Vector<scalar_t, lno_t, gno_t, node_t> >
759 typedef Xpetra::Vector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
760 typedef Xpetra::TpetraVector<scalar_t, lno_t, gno_t, node_t> xt_vector_t;
761 typedef Tpetra::Vector<scalar_t, lno_t, gno_t, node_t> t_vector_t;
763 static inline RCP<x_vector_t>
convertToXpetra(
const RCP<x_vector_t> &a)
768 static RCP<x_vector_t>
doMigration(
const x_vector_t &from,
769 size_t numLocalRows,
const gno_t *myNewRows)
771 Xpetra::UnderlyingLib lib = from.getMap()->lib();
773 if (lib == Xpetra::UseEpetra){
774 throw std::logic_error(
"compiler should have used specialization");
777 const xt_vector_t *xtv =
dynamic_cast<const xt_vector_t *
>(&from);
778 RCP<const t_vector_t> tv = xtv->getTpetra_Vector();
781 *tv, numLocalRows, myNewRows);
792 template <
typename node_t>
793 struct XpetraTraits<Xpetra::
Vector<double, int, int, node_t> >
795 typedef double scalar_t;
798 typedef Xpetra::Vector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
799 typedef Xpetra::TpetraVector<scalar_t, lno_t, gno_t, node_t> xt_vector_t;
800 typedef Tpetra::Vector<scalar_t, lno_t, gno_t, node_t> t_vector_t;
802 static inline RCP<x_vector_t>
convertToXpetra(
const RCP<x_vector_t> &a)
807 static RCP<x_vector_t>
doMigration(
const x_vector_t &from,
808 size_t numLocalRows,
const gno_t *myNewRows)
810 Xpetra::UnderlyingLib lib = from.getMap()->lib();
812 if (lib == Xpetra::UseEpetra){
813 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA)
814 typedef Epetra_Vector e_vector_t;
815 typedef Xpetra::EpetraVectorT<gno_t,node_t> xe_vector_t;
817 const xe_vector_t *xev =
dynamic_cast<const xe_vector_t *
>(&from);
818 RCP<const e_vector_t> ev = rcp(xev->getEpetra_Vector());
821 *ev, numLocalRows, myNewRows);
827 throw std::runtime_error(
"Xpetra with Epetra requested, but "
828 "Trilinos is not built with Epetra");
832 const xt_vector_t *xtv =
dynamic_cast<const xt_vector_t *
>(&from);
833 RCP<t_vector_t> tv = xtv->getTpetra_Vector();
836 *tv, numLocalRows, myNewRows);
847 template <
typename scalar_t,
851 struct XpetraTraits<Tpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> >
853 typedef Tpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> t_vector_t;
854 typedef Xpetra::TpetraMultiVector<scalar_t, lno_t, gno_t, node_t> xt_vector_t;
855 typedef Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
857 static inline RCP<x_vector_t>
convertToXpetra(
const RCP<t_vector_t> &a)
859 return rcp(
new xt_vector_t(a));
862 static RCP<t_vector_t>
doMigration(
const t_vector_t &from,
863 size_t numLocalElts,
const gno_t *myNewElts)
865 typedef Tpetra::Map<lno_t, gno_t, node_t>
map_t;
868 const RCP<const map_t> &smap = from.getMap();
869 gno_t numGlobalElts = smap->getGlobalNumElements();
870 gno_t base = smap->getMinAllGlobalIndex();
873 ArrayView<const gno_t> eltList(myNewElts, numLocalElts);
874 const RCP<const Teuchos::Comm<int> > comm = from.getMap()->getComm();
875 RCP<const map_t> tmap = rcp(
new map_t(numGlobalElts, eltList, base, comm));
878 Tpetra::Import<lno_t, gno_t, node_t> importer(smap, tmap);
881 RCP<t_vector_t> MV = rcp(
882 new t_vector_t(tmap, from.getNumVectors(),
true));
883 MV->doImport(from, importer, Tpetra::INSERT);
890 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA)
893 struct XpetraTraits<Epetra_MultiVector>
898 typedef InputTraits<Epetra_MultiVector>::scalar_t scalar_t;
899 typedef Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> x_mvector_t;
902 const RCP<Epetra_MultiVector> &a)
904 RCP<Xpetra::EpetraMultiVectorT<gno_t, node_t> > xemv;
906 xemv = rcp(
new Xpetra::EpetraMultiVectorT<gno_t,node_t>(a));
908 catch (std::exception &e) {
909 if (std::is_same<node_t, Xpetra::EpetraNode>::value)
910 throw std::runtime_error(std::string(
"Cannot convert from "
911 "Epetra_MultiVector to "
912 "Xpetra::EpetraMultiVectorT\n")
915 throw std::runtime_error(std::string(
"Cannot convert from "
916 "Epetra_MultiVector to "
917 "Xpetra::EpetraMultiVectorT\n"
918 "Use node_t that is supported by "
919 "Xpetra with Epetra classes\n")
922 return rcp_implicit_cast<x_mvector_t>(xemv);
925 static RCP<Epetra_MultiVector>
doMigration(
const Epetra_MultiVector &from,
926 size_t numLocalElts,
const gno_t *myNewElts)
929 const Epetra_BlockMap &smap = from.Map();
930 gno_t numGlobalElts = smap.NumGlobalElements();
931 int base = smap.MinAllGID();
934 const Epetra_Comm &comm = from.Comm();
935 const Epetra_BlockMap tmap(numGlobalElts, numLocalElts, myNewElts,
939 Epetra_Import importer(tmap, smap);
942 RCP<Epetra_MultiVector> MV = rcp(
943 new Epetra_MultiVector(tmap, from.NumVectors(),
true));
944 Epetra_CombineMode c = Insert;
945 MV->Import(from, importer, c);
954 template <
typename scalar_t,
958 struct XpetraTraits<Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> >
960 typedef Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> x_mvector_t;
961 typedef Xpetra::TpetraMultiVector<scalar_t, lno_t, gno_t, node_t> xt_mvector_t;
962 typedef Tpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> t_mvector_t;
964 static inline RCP<x_mvector_t>
convertToXpetra(
const RCP<x_mvector_t> &a)
969 static RCP<x_mvector_t>
doMigration(
const x_mvector_t &from,
970 size_t numLocalRows,
const gno_t *myNewRows)
972 Xpetra::UnderlyingLib lib = from.getMap()->lib();
974 if (lib == Xpetra::UseEpetra){
975 throw std::logic_error(
"compiler should have used specialization");
978 const xt_mvector_t *xtv =
dynamic_cast<const xt_mvector_t *
>(&from);
979 RCP<t_mvector_t> tv = xtv->getTpetra_MultiVector();
982 *tv, numLocalRows, myNewRows);
993 template <
typename node_t>
994 struct XpetraTraits<Xpetra::MultiVector<double, int, int, node_t> >
996 typedef double scalar_t;
999 typedef Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> x_mvector_t;
1000 typedef Xpetra::TpetraMultiVector<scalar_t, lno_t, gno_t, node_t> xt_mvector_t;
1001 typedef Tpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> t_mvector_t;
1003 static inline RCP<x_mvector_t>
convertToXpetra(
const RCP<x_mvector_t> &a)
1008 static RCP<x_mvector_t>
doMigration(
const x_mvector_t &from,
1009 size_t numLocalRows,
const gno_t *myNewRows)
1011 Xpetra::UnderlyingLib lib = from.getMap()->lib();
1013 if (lib == Xpetra::UseEpetra){
1014 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA)
1015 typedef Epetra_MultiVector e_mvector_t;
1016 typedef Xpetra::EpetraMultiVectorT<gno_t,node_t> xe_mvector_t;
1018 const xe_mvector_t *xev =
dynamic_cast<const xe_mvector_t *
>(&from);
1019 RCP<e_mvector_t> ev = xev->getEpetra_MultiVector();
1022 *ev, numLocalRows, myNewRows);
1028 throw std::runtime_error(
"Xpetra with Epetra requested, but "
1029 "Trilinos is not built with Epetra");
1033 const xt_mvector_t *xtv =
dynamic_cast<const xt_mvector_t *
>(&from);
1034 RCP<t_mvector_t> tv = xtv->getTpetra_MultiVector();
1037 *tv, numLocalRows, myNewRows);
1046 #endif // DOXYGEN_SHOULD_SKIP_THIS
1050 #endif // _ZOLTAN2_XPETRATRAITS_HPP_
Defines the traits required for Tpetra, Eptra and Xpetra objects.
default_gno_t gno_t
The objects global ordinal data type.
map_t::global_ordinal_type gno_t
static RCP< User > doMigration(const User &from, size_t numLocalRows, const gno_t *myNewRows)
Migrate the object Given a user object and a new row distribution, create and return a new user objec...
::Tpetra::Details::DefaultTypes::global_ordinal_type default_gno_t
::Tpetra::Details::DefaultTypes::local_ordinal_type default_lno_t
Tpetra::Vector< z2TestScalar, z2TestLO, z2TestGO > Vector
static RCP< User > convertToXpetra(const RCP< User > &a)
Convert the object to its Xpetra wrapped version.
map_t::local_ordinal_type lno_t
default_lno_t lno_t
The objects local ordinal data type.
Gathering definitions used in software development.