50 #ifndef _ZOLTAN2_XPETRATRAITS_HPP_ 51 #define _ZOLTAN2_XPETRATRAITS_HPP_ 56 #include <Xpetra_TpetraCrsMatrix.hpp> 57 #include <Xpetra_TpetraRowMatrix.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>
136 typedef typename Xpetra::TpetraCrsMatrix<scalar_t,lno_t,gno_t,node_t>
138 typedef typename Tpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t>
143 return rcp(
new xtmatrix_t(a));
146 static RCP<tmatrix_t>
doMigration(
const tmatrix_t &from,
147 size_t numLocalRows,
const gno_t *myNewRows)
149 typedef Tpetra::Map<lno_t, gno_t, node_t> map_t;
153 const RCP<const map_t> &smap = from.getRowMap();
154 int oldNumElts = smap->getNodeNumElements();
155 gno_t numGlobalRows = smap->getGlobalNumElements();
158 ArrayView<const gno_t> rowList(myNewRows, numLocalRows);
159 const RCP<const Teuchos::Comm<int> > &comm = from.getComm();
160 RCP<const map_t> tmap = rcp(
161 new map_t(numGlobalRows, rowList, base, comm));
162 int newNumElts = numLocalRows;
165 Tpetra::Import<lno_t, gno_t, node_t> importer(smap, tmap);
168 typedef Tpetra::Vector<scalar_t, lno_t, gno_t, node_t> vector_t;
169 vector_t numOld(smap);
170 vector_t numNew(tmap);
171 for (
int lid=0; lid < oldNumElts; lid++){
172 numOld.replaceGlobalValue(smap->getGlobalElement(lid),
173 scalar_t(from.getNumEntriesInLocalRow(lid)));
175 numNew.doImport(numOld, importer, Tpetra::INSERT);
178 ArrayRCP<size_t> nnz(newNumElts);
180 ArrayRCP<scalar_t> ptr = numNew.getDataNonConst(0);
181 for (
int lid=0; lid < newNumElts; lid++){
182 nnz[lid] =
static_cast<size_t>(ptr[lid]);
187 RCP<tmatrix_t> M = rcp(
new tmatrix_t(tmap, nnz, Tpetra::StaticProfile));
188 M->doImport(from, importer, Tpetra::INSERT);
196 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA) 206 static inline RCP<Xpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t> >
209 return rcp(
new Xpetra::EpetraCrsMatrixT<gno_t, node_t>(a));
213 static RCP<Epetra_CrsMatrix>
doMigration(
const Epetra_CrsMatrix &from,
214 size_t numLocalRows,
const gno_t *myNewRows)
219 const Epetra_Map &smap = from.RowMap();
220 int oldNumElts = smap.NumMyElements();
221 gno_t numGlobalRows = smap.NumGlobalElements();
224 const Epetra_Comm &comm = from.Comm();
225 Epetra_Map tmap(numGlobalRows, numLocalRows, myNewRows, base, comm);
226 int newNumElts = numLocalRows;
229 Epetra_Import importer(tmap, smap);
232 Epetra_Vector numOld(smap);
233 Epetra_Vector numNew(tmap);
235 for (
int lid=0; lid < oldNumElts; lid++){
236 numOld[lid] = from.NumMyEntries(lid);
238 numNew.Import(numOld, importer, Insert);
240 Array<int> nnz(newNumElts);
241 for (
int lid=0; lid < newNumElts; lid++){
242 nnz[lid] =
static_cast<int>(numNew[lid]);
246 RCP<Epetra_CrsMatrix> M = rcp(
new Epetra_CrsMatrix(::Copy, tmap, nnz.getRawPtr(),
true));
247 M->Import(from, importer, Insert);
259 template <
typename scalar_t,
263 struct XpetraTraits<Xpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> >
265 typedef Xpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> x_matrix_t;
266 typedef Xpetra::TpetraCrsMatrix<scalar_t, lno_t, gno_t, node_t> xt_matrix_t;
267 typedef Tpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t> t_matrix_t;
269 static inline RCP<x_matrix_t>
convertToXpetra(
const RCP<x_matrix_t > &a)
274 static RCP<x_matrix_t>
doMigration(
const x_matrix_t &from,
275 size_t numLocalRows,
const gno_t *myNewRows)
277 Xpetra::UnderlyingLib lib = from.getRowMap()->lib();
279 if (lib == Xpetra::UseEpetra){
280 throw std::logic_error(
"compiler should have used specialization");
283 const xt_matrix_t *xtm =
dynamic_cast<const xt_matrix_t *
>(&from);
284 RCP<const t_matrix_t> tm = xtm->getTpetra_CrsMatrix();
287 *tm, numLocalRows, myNewRows);
299 template <
typename node_t>
300 struct XpetraTraits<Xpetra::CrsMatrix<double, int, int, node_t> >
302 typedef double scalar_t;
305 typedef Xpetra::CrsMatrix<scalar_t, lno_t, gno_t, node_t> x_matrix_t;
306 typedef Xpetra::TpetraCrsMatrix<scalar_t, lno_t, gno_t, node_t> xt_matrix_t;
307 typedef Tpetra::CrsMatrix<scalar_t,lno_t,gno_t,node_t> t_matrix_t;
309 static inline RCP<x_matrix_t>
convertToXpetra(
const RCP<x_matrix_t > &a)
314 static RCP<x_matrix_t>
doMigration(
const x_matrix_t &from,
315 size_t numLocalRows,
const gno_t *myNewRows)
317 Xpetra::UnderlyingLib lib = from.getRowMap()->lib();
319 if (lib == Xpetra::UseEpetra){
320 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA) 321 typedef Epetra_CrsMatrix e_matrix_t;
322 typedef Xpetra::EpetraCrsMatrixT<gno_t,node_t> xe_matrix_t;
324 const xe_matrix_t *xem =
dynamic_cast<const xe_matrix_t *
>(&from);
325 RCP<const e_matrix_t> em = xem->getEpetra_CrsMatrix();
328 *em, numLocalRows, myNewRows);
334 throw std::runtime_error(
"Xpetra with Epetra requested, but " 335 "Trilinos is not built with Epetra");
339 const xt_matrix_t *xtm =
dynamic_cast<const xt_matrix_t *
>(&from);
340 RCP<const t_matrix_t> tm = xtm->getTpetra_CrsMatrix();
343 *tm, numLocalRows, myNewRows);
355 template <
typename lno_t,
358 struct XpetraTraits<Tpetra::CrsGraph<lno_t, gno_t, node_t> >
360 typedef typename Xpetra::CrsGraph<lno_t, gno_t, node_t>
xgraph_t;
361 typedef typename Xpetra::TpetraCrsGraph<lno_t, gno_t, node_t> xtgraph_t;
362 typedef typename Tpetra::CrsGraph<lno_t, gno_t, node_t>
tgraph_t;
366 return rcp(
new xtgraph_t(a));
369 static RCP<tgraph_t>
doMigration(
const tgraph_t &from,
370 size_t numLocalRows,
const gno_t *myNewRows)
372 typedef Tpetra::Map<lno_t, gno_t, node_t> map_t;
376 const RCP<const map_t> &smap = from.getRowMap();
377 int oldNumElts = smap->getNodeNumElements();
378 gno_t numGlobalRows = smap->getGlobalNumElements();
381 ArrayView<const gno_t> rowList(myNewRows, numLocalRows);
382 const RCP<const Teuchos::Comm<int> > &comm = from.getComm();
383 RCP<const map_t> tmap = rcp(
384 new map_t(numGlobalRows, rowList, base, comm));
387 Tpetra::Import<lno_t, gno_t, node_t> importer(smap, tmap);
390 typedef Tpetra::Vector<gno_t, lno_t, gno_t, node_t> vector_t;
391 vector_t numOld(smap);
392 vector_t numNew(tmap);
393 for (
int lid=0; lid < oldNumElts; lid++){
394 numOld.replaceGlobalValue(smap->getGlobalElement(lid),
395 from.getNumEntriesInLocalRow(lid));
397 numNew.doImport(numOld, importer, Tpetra::INSERT);
399 size_t numElts = tmap->getNodeNumElements();
400 ArrayRCP<const gno_t> nnz;
402 nnz = numNew.getData(0);
404 ArrayRCP<const size_t> nnz_size_t;
406 if (numElts &&
sizeof(gno_t) !=
sizeof(
size_t)){
407 size_t *
vals =
new size_t [numElts];
408 nnz_size_t = arcp(vals, 0, numElts,
true);
409 for (
size_t i=0; i < numElts; i++){
410 vals[i] =
static_cast<size_t>(nnz[i]);
414 nnz_size_t = arcp_reinterpret_cast<
const size_t>(nnz);
418 RCP<tgraph_t> G = rcp(
new tgraph_t(tmap, nnz_size_t, Tpetra::StaticProfile));
420 G->doImport(from, importer, Tpetra::INSERT);
429 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA) 437 static inline RCP<Xpetra::CrsGraph<lno_t,gno_t,node_t> >
440 return rcp(
new Xpetra::EpetraCrsGraphT<gno_t,node_t>(a));
443 static RCP<Epetra_CrsGraph>
doMigration(
const Epetra_CrsGraph &from,
444 size_t numLocalRows,
const gno_t *myNewRows)
449 const Epetra_BlockMap &smap = from.RowMap();
450 gno_t numGlobalRows = smap.NumGlobalElements();
451 lno_t oldNumElts = smap.NumMyElements();
454 const Epetra_Comm &comm = from.Comm();
455 Epetra_BlockMap tmap(numGlobalRows, numLocalRows,
456 myNewRows, 1, base, comm);
457 lno_t newNumElts = tmap.NumMyElements();
460 Epetra_Import importer(tmap, smap);
463 Epetra_Vector numOld(smap);
464 Epetra_Vector numNew(tmap);
466 for (
int lid=0; lid < oldNumElts; lid++){
467 numOld[lid] = from.NumMyIndices(lid);
469 numNew.Import(numOld, importer, Insert);
471 Array<int> nnz(newNumElts);
472 for (
int lid=0; lid < newNumElts; lid++){
473 nnz[lid] =
static_cast<int>(numNew[lid]);
477 RCP<Epetra_CrsGraph> G = rcp(
new Epetra_CrsGraph(::Copy, tmap, nnz.getRawPtr(),
true));
478 G->Import(from, importer, Insert);
491 template <
typename lno_t,
494 struct XpetraTraits<Xpetra::CrsGraph<lno_t, gno_t, node_t> >
496 typedef Xpetra::CrsGraph<lno_t, gno_t, node_t> x_graph_t;
497 typedef Xpetra::TpetraCrsGraph<lno_t, gno_t, node_t> xt_graph_t;
498 typedef Tpetra::CrsGraph<lno_t,gno_t,node_t> t_graph_t;
505 static RCP<x_graph_t>
doMigration(
const x_graph_t &from,
506 size_t numLocalRows,
const gno_t *myNewRows)
508 Xpetra::UnderlyingLib lib = from.getRowMap()->lib();
510 if (lib == Xpetra::UseEpetra){
511 throw std::logic_error(
"compiler should have used specialization");
514 const xt_graph_t *xtg =
dynamic_cast<const xt_graph_t *
>(&from);
515 RCP<const t_graph_t> tg = xtg->getTpetra_CrsGraph();
518 *tg, numLocalRows, myNewRows);
530 template <
typename scalar_t,
534 struct XpetraTraits<Xpetra::RowMatrix<scalar_t, lno_t, gno_t, node_t> >
536 typedef Xpetra::RowMatrix<scalar_t, lno_t, gno_t, node_t> x_matrix_t;
537 typedef Xpetra::TpetraRowMatrix<scalar_t, lno_t, gno_t, node_t> xt_matrix_t;
538 typedef Tpetra::RowMatrix<scalar_t,lno_t,gno_t,node_t> t_matrix_t;
540 static inline RCP<x_matrix_t>
convertToXpetra(
const RCP<x_matrix_t > &a)
545 static RCP<x_matrix_t>
doMigration(
const x_matrix_t &from,
546 size_t numLocalRows,
const gno_t *myNewRows)
548 Xpetra::UnderlyingLib lib = from.getRowMap()->lib();
550 if (lib == Xpetra::UseEpetra){
551 throw std::logic_error(
"compiler should have used specialization");
554 const xt_matrix_t *xtm =
dynamic_cast<const xt_matrix_t *
>(&from);
555 RCP<const t_matrix_t> tm = xtm->getTpetra_CrsMatrix();
558 *tm, numLocalRows, myNewRows);
569 template <
typename node_t>
570 struct XpetraTraits<Xpetra::CrsGraph<int, int, node_t> >
574 typedef Xpetra::CrsGraph<lno_t, gno_t, node_t> x_graph_t;
575 typedef Xpetra::TpetraCrsGraph<lno_t, gno_t, node_t> xt_graph_t;
576 typedef Tpetra::CrsGraph<lno_t,gno_t,node_t> t_graph_t;
583 static RCP<x_graph_t>
doMigration(
const x_graph_t &from,
584 size_t numLocalRows,
const gno_t *myNewRows)
586 Xpetra::UnderlyingLib lib = from.getRowMap()->lib();
588 if (lib == Xpetra::UseEpetra){
589 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA) 590 typedef Xpetra::EpetraCrsGraphT<gno_t,node_t> xe_graph_t;
591 typedef Epetra_CrsGraph e_graph_t;
593 const xe_graph_t *xeg =
dynamic_cast<const xe_graph_t *
>(&from);
594 RCP<const e_graph_t> eg = xeg->getEpetra_CrsGraph();
597 *eg, numLocalRows, myNewRows);
603 throw std::runtime_error(
"Xpetra with Epetra requested, but " 604 "Trilinos is not built with Epetra");
608 const xt_graph_t *xtg =
dynamic_cast<const xt_graph_t *
>(&from);
609 RCP<const t_graph_t> tg = xtg->getTpetra_CrsGraph();
612 *tg, numLocalRows, myNewRows);
623 template <
typename scalar_t,
629 typedef Tpetra::Vector<scalar_t, lno_t, gno_t, node_t> t_vector_t;
630 typedef Xpetra::TpetraVector<scalar_t, lno_t, gno_t, node_t> xt_vector_t;
631 typedef Xpetra::Vector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
633 static inline RCP<x_vector_t>
convertToXpetra(
const RCP<t_vector_t> &a)
635 return rcp(
new xt_vector_t(a));
638 static RCP<t_vector_t>
doMigration(
const t_vector_t &from,
639 size_t numLocalElts,
const gno_t *myNewElts)
642 typedef Tpetra::Map<lno_t, gno_t, node_t> map_t;
645 const RCP<const map_t> &smap = from.getMap();
646 gno_t numGlobalElts = smap->getGlobalNumElements();
649 ArrayView<const gno_t> eltList(myNewElts, numLocalElts);
650 const RCP<const Teuchos::Comm<int> > comm = from.getMap()->getComm();
651 RCP<const map_t> tmap = rcp(
652 new map_t(numGlobalElts, eltList, base, comm));
655 Tpetra::Import<lno_t, gno_t, node_t> importer(smap, tmap);
659 Tpetra::createVector<scalar_t,lno_t,gno_t,node_t>(tmap);
660 V->doImport(from, importer, Tpetra::INSERT);
667 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA) 677 typedef Xpetra::Vector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
679 static inline RCP<x_vector_t>
convertToXpetra(
const RCP<Epetra_Vector> &a)
681 RCP<Xpetra::EpetraVectorT<gno_t,node_t> > xev = rcp(
new Xpetra::EpetraVectorT<gno_t,node_t>(a));
682 return rcp_implicit_cast<x_vector_t>(xev);
685 static RCP<Epetra_Vector>
doMigration(
const Epetra_Vector &from,
686 size_t numLocalElts,
const gno_t *myNewElts)
690 const Epetra_BlockMap &smap = from.Map();
691 gno_t numGlobalElts = smap.NumGlobalElements();
694 const Epetra_Comm &comm = from.Comm();
695 const Epetra_BlockMap tmap(numGlobalElts, numLocalElts, myNewElts,
699 Epetra_Import importer(tmap, smap);
702 RCP<Epetra_Vector> V = rcp(
new Epetra_Vector(tmap,
true));
703 Epetra_CombineMode c = Insert;
704 V->Import(from, importer, c);
713 template <
typename scalar_t,
719 typedef Xpetra::Vector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
720 typedef Xpetra::TpetraVector<scalar_t, lno_t, gno_t, node_t> xt_vector_t;
721 typedef Tpetra::Vector<scalar_t, lno_t, gno_t, node_t> t_vector_t;
723 static inline RCP<x_vector_t>
convertToXpetra(
const RCP<x_vector_t> &a)
728 static RCP<x_vector_t>
doMigration(
const x_vector_t &from,
729 size_t numLocalRows,
const gno_t *myNewRows)
731 Xpetra::UnderlyingLib lib = from.getMap()->lib();
733 if (lib == Xpetra::UseEpetra){
734 throw std::logic_error(
"compiler should have used specialization");
737 const xt_vector_t *xtv =
dynamic_cast<const xt_vector_t *
>(&from);
738 RCP<const t_vector_t> tv = xtv->getTpetra_Vector();
741 *tv, numLocalRows, myNewRows);
752 template <
typename node_t>
755 typedef double scalar_t;
758 typedef Xpetra::Vector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
759 typedef Xpetra::TpetraVector<scalar_t, lno_t, gno_t, node_t> xt_vector_t;
760 typedef Tpetra::Vector<scalar_t, lno_t, gno_t, node_t> t_vector_t;
762 static inline RCP<x_vector_t>
convertToXpetra(
const RCP<x_vector_t> &a)
767 static RCP<x_vector_t>
doMigration(
const x_vector_t &from,
768 size_t numLocalRows,
const gno_t *myNewRows)
770 Xpetra::UnderlyingLib lib = from.getMap()->lib();
772 if (lib == Xpetra::UseEpetra){
773 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA) 774 typedef Epetra_Vector e_vector_t;
775 typedef Xpetra::EpetraVectorT<gno_t,node_t> xe_vector_t;
777 const xe_vector_t *xev =
dynamic_cast<const xe_vector_t *
>(&from);
778 RCP<const e_vector_t> ev = rcp(xev->getEpetra_Vector());
781 *ev, numLocalRows, myNewRows);
787 throw std::runtime_error(
"Xpetra with Epetra requested, but " 788 "Trilinos is not built with Epetra");
792 const xt_vector_t *xtv =
dynamic_cast<const xt_vector_t *
>(&from);
793 RCP<t_vector_t> tv = xtv->getTpetra_Vector();
796 *tv, numLocalRows, myNewRows);
807 template <
typename scalar_t,
811 struct XpetraTraits<Tpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> >
813 typedef Tpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> t_vector_t;
814 typedef Xpetra::TpetraMultiVector<scalar_t, lno_t, gno_t, node_t> xt_vector_t;
815 typedef Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> x_vector_t;
817 static inline RCP<x_vector_t>
convertToXpetra(
const RCP<t_vector_t> &a)
819 return rcp(
new xt_vector_t(a));
822 static RCP<t_vector_t>
doMigration(
const t_vector_t &from,
823 size_t numLocalElts,
const gno_t *myNewElts)
825 typedef Tpetra::Map<lno_t, gno_t, node_t> map_t;
829 const RCP<const map_t> &smap = from.getMap();
830 gno_t numGlobalElts = smap->getGlobalNumElements();
833 ArrayView<const gno_t> eltList(myNewElts, numLocalElts);
834 const RCP<const Teuchos::Comm<int> > comm = from.getMap()->getComm();
835 RCP<const map_t> tmap = rcp(
836 new map_t(numGlobalElts, eltList, base, comm));
839 Tpetra::Import<lno_t, gno_t, node_t> importer(smap, tmap);
842 RCP<t_vector_t> MV = rcp(
843 new t_vector_t(tmap, from.getNumVectors(),
true));
844 MV->doImport(from, importer, Tpetra::INSERT);
851 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA) 860 typedef Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> x_mvector_t;
863 const RCP<Epetra_MultiVector> &a)
865 RCP<Xpetra::EpetraMultiVectorT<gno_t,node_t> > xemv = rcp(
new Xpetra::EpetraMultiVectorT<gno_t,node_t>(a));
866 return rcp_implicit_cast<x_mvector_t>(xemv);
869 static RCP<Epetra_MultiVector>
doMigration(
const Epetra_MultiVector &from,
870 size_t numLocalElts,
const gno_t *myNewElts)
874 const Epetra_BlockMap &smap = from.Map();
875 gno_t numGlobalElts = smap.NumGlobalElements();
878 const Epetra_Comm &comm = from.Comm();
879 const Epetra_BlockMap tmap(numGlobalElts, numLocalElts, myNewElts,
883 Epetra_Import importer(tmap, smap);
886 RCP<Epetra_MultiVector> MV = rcp(
887 new Epetra_MultiVector(tmap, from.NumVectors(),
true));
888 Epetra_CombineMode c = Insert;
889 MV->Import(from, importer, c);
898 template <
typename scalar_t,
902 struct XpetraTraits<Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> >
904 typedef Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> x_mvector_t;
905 typedef Xpetra::TpetraMultiVector<scalar_t, lno_t, gno_t, node_t> xt_mvector_t;
906 typedef Tpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> t_mvector_t;
908 static inline RCP<x_mvector_t>
convertToXpetra(
const RCP<x_mvector_t> &a)
913 static RCP<x_mvector_t>
doMigration(
const x_mvector_t &from,
914 size_t numLocalRows,
const gno_t *myNewRows)
916 Xpetra::UnderlyingLib lib = from.getMap()->lib();
918 if (lib == Xpetra::UseEpetra){
919 throw std::logic_error(
"compiler should have used specialization");
922 const xt_mvector_t *xtv =
dynamic_cast<const xt_mvector_t *
>(&from);
923 RCP<t_mvector_t> tv = xtv->getTpetra_MultiVector();
926 *tv, numLocalRows, myNewRows);
937 template <
typename node_t>
938 struct XpetraTraits<Xpetra::MultiVector<double, int, int, node_t> >
940 typedef double scalar_t;
943 typedef Xpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> x_mvector_t;
944 typedef Xpetra::TpetraMultiVector<scalar_t, lno_t, gno_t, node_t> xt_mvector_t;
945 typedef Tpetra::MultiVector<scalar_t, lno_t, gno_t, node_t> t_mvector_t;
947 static inline RCP<x_mvector_t>
convertToXpetra(
const RCP<x_mvector_t> &a)
952 static RCP<x_mvector_t>
doMigration(
const x_mvector_t &from,
953 size_t numLocalRows,
const gno_t *myNewRows)
955 Xpetra::UnderlyingLib lib = from.getMap()->lib();
957 if (lib == Xpetra::UseEpetra){
958 #if defined(HAVE_ZOLTAN2_EPETRA) && defined(HAVE_XPETRA_EPETRA) 959 typedef Epetra_MultiVector e_mvector_t;
960 typedef Xpetra::EpetraMultiVectorT<gno_t,node_t> xe_mvector_t;
962 const xe_mvector_t *xev =
dynamic_cast<const xe_mvector_t *
>(&from);
963 RCP<e_mvector_t> ev = xev->getEpetra_MultiVector();
966 *ev, numLocalRows, myNewRows);
972 throw std::runtime_error(
"Xpetra with Epetra requested, but " 973 "Trilinos is not built with Epetra");
977 const xt_mvector_t *xtv =
dynamic_cast<const xt_mvector_t *
>(&from);
978 RCP<t_mvector_t> tv = xtv->getTpetra_MultiVector();
981 *tv, numLocalRows, myNewRows);
990 #endif // DOXYGEN_SHOULD_SKIP_THIS 994 #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.
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...
static RCP< User > convertToXpetra(const RCP< User > &a)
Convert the object to its Xpetra wrapped version.
default_lno_t lno_t
The objects local ordinal data type.
Tpetra::Vector< z2TestScalar, z2TestLO, z2TestGO > Vector
Gathering definitions used in software development.