46#ifndef MUELU_UTILITIES_DEF_HPP
47#define MUELU_UTILITIES_DEF_HPP
49#include <Teuchos_DefaultComm.hpp>
50#include <Teuchos_ParameterList.hpp>
54#ifdef HAVE_MUELU_EPETRA
56# include "Epetra_MpiComm.h"
60#if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT)
67#include <Xpetra_EpetraUtils.hpp>
68#include <Xpetra_EpetraMultiVector.hpp>
72#ifdef HAVE_MUELU_TPETRA
73#include <MatrixMarket_Tpetra.hpp>
74#include <Tpetra_RowMatrixTransposer.hpp>
75#include <TpetraExt_MatrixMatrix.hpp>
76#include <Xpetra_TpetraMultiVector.hpp>
77#include <Xpetra_TpetraOperator.hpp>
78#include <Xpetra_TpetraCrsMatrix.hpp>
79#include <Xpetra_TpetraBlockCrsMatrix.hpp>
82#ifdef HAVE_MUELU_EPETRA
83#include <Xpetra_EpetraMap.hpp>
86#include <Xpetra_BlockedCrsMatrix.hpp>
88#include <Xpetra_IO.hpp>
89#include <Xpetra_Import.hpp>
90#include <Xpetra_ImportFactory.hpp>
91#include <Xpetra_Map.hpp>
92#include <Xpetra_MapFactory.hpp>
93#include <Xpetra_Matrix.hpp>
94#include <Xpetra_MatrixFactory.hpp>
95#include <Xpetra_MultiVector.hpp>
96#include <Xpetra_MultiVectorFactory.hpp>
97#include <Xpetra_Operator.hpp>
98#include <Xpetra_Vector.hpp>
99#include <Xpetra_VectorFactory.hpp>
101#include <Xpetra_MatrixMatrix.hpp>
104#if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_ML)
105#include <ml_operator.h>
106#include <ml_epetra_utils.h>
111#ifdef HAVE_MUELU_EPETRA
116#ifdef HAVE_MUELU_EPETRA
117 template<
typename SC,
typename LO,
typename GO,
typename NO>
119 return Xpetra::Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap<SC,LO,GO,NO>(epAB);
123#ifdef HAVE_MUELU_EPETRA
124 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
127 if (
tmpVec == Teuchos::null)
128 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
129 return tmpVec->getEpetra_MultiVector();
132 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
135 if (
tmpVec == Teuchos::null)
136 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
137 return tmpVec->getEpetra_MultiVector();
140 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
142 const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> &
tmpVec =
dynamic_cast<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> &
>(
vec);
143 return *(
tmpVec.getEpetra_MultiVector());
146 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
148 const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> &
tmpVec =
dynamic_cast<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> &
>(
vec);
149 return *(
tmpVec.getEpetra_MultiVector());
152 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
155 if (
crsOp == Teuchos::null)
159 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
163 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
166 if (
crsOp == Teuchos::null)
170 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
174 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
177 const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
crsOp =
dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(Op);
179 const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>&
tmp_ECrsMtx =
dynamic_cast<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>&
>(*
crsOp.getCrsMatrix());
181 }
catch (std::bad_cast&) {
182 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
184 }
catch (std::bad_cast&) {
189 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
192 Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
crsOp =
dynamic_cast<Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(Op);
194 Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>&
tmp_ECrsMtx =
dynamic_cast<Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>&
>(*
crsOp.getCrsMatrix());
196 }
catch (std::bad_cast&) {
197 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
199 }
catch (std::bad_cast&) {
204 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
207 if (
xeMap == Teuchos::null)
208 throw Exceptions::BadCast(
"Utilities::Map2EpetraMap : Cast from Xpetra::Map to Xpetra::EpetraMap failed");
209 return xeMap->getEpetra_Map();
213#ifdef HAVE_MUELU_TPETRA
214 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
215 RCP<const Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
218 if (
tmpVec == Teuchos::null)
219 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
220 return tmpVec->getTpetra_MultiVector();
223 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
226 if (
tmpVec == Teuchos::null)
227 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
228 return tmpVec->getTpetra_MultiVector();
231 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
233 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
tmpVec =
dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(
vec);
234 return *(
tmpVec.getTpetra_MultiVector());
237 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
239 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
tmpVec =
dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(
vec);
240 return tmpVec.getTpetra_MultiVector();
243 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
244 const Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>&
246 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
tmpVec =
dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(
vec);
247 return *(
tmpVec.getTpetra_MultiVector());
250 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
254 if (
crsOp == Teuchos::null)
258 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
262 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
265 if (
crsOp == Teuchos::null)
269 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
273 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
276 const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
crsOp =
dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(Op);
278 const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
tmp_ECrsMtx =
dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(*
crsOp.getCrsMatrix());
280 }
catch (std::bad_cast&) {
281 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
283 }
catch (std::bad_cast&) {
288 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
291 Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
crsOp =
dynamic_cast<Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(Op);
293 Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
tmp_ECrsMtx =
dynamic_cast<Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(*
crsOp.getCrsMatrix());
295 }
catch (std::bad_cast&) {
296 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
298 }
catch (std::bad_cast&) {
304 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
306 using XCrsMatrixWrap = Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
309 if (
crsOp == Teuchos::null)
313 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
317 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
319 using XCrsMatrixWrap = Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
321 if (
crsOp == Teuchos::null)
325 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
326 return tmp_ECrsMtx->getTpetra_BlockCrsMatrixNonConst();
329 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
332 using XCrsMatrixWrap = Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
335 const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
tmp_ECrsMtx =
dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(*
crsOp.getCrsMatrix());
337 }
catch (std::bad_cast&) {
338 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
340 }
catch (std::bad_cast&) {
345 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
348 using XCrsMatrixWrap = Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
351 Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
tmp_ECrsMtx =
dynamic_cast<Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(*
crsOp.getCrsMatrix());
352 return *
tmp_ECrsMtx.getTpetra_BlockCrsMatrixNonConst();
353 }
catch (std::bad_cast&) {
354 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
356 }
catch (std::bad_cast&) {
362 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
365 if (!
mat.is_null()) {
367 if (
crsOp == Teuchos::null)
374 return tmp_Crs->getTpetra_CrsMatrixNonConst();
379 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
380 return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
390 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
393 if (!
mat.is_null()) {
395 if (
crsOp == Teuchos::null)
402 return tmp_Crs->getTpetra_CrsMatrixNonConst();
407 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
408 return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
419 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
423 throw Exceptions::BadCast(
"Utilities::Map2TpetraMap : Cast from Xpetra::Map to Xpetra::TpetraMap failed");
428 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
433 Scalar one = Teuchos::ScalarTraits<Scalar>::one();
443 switch (Op.getRowMap()->lib()) {
444 case Xpetra::UseTpetra:
448 case Xpetra::UseEpetra:
457 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
462 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
467#ifdef HAVE_MUELU_TPETRA
469 Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>&
tpOp = Op2NonConstTpetraCrs(Op);
480 if (
tpOp.isFillComplete())
483 if (Op.isLocallyIndexed() ==
true) {
484 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::local_inds_host_view_type
cols;
485 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::values_host_view_type
vals;
487 for (
size_t i = 0; i < rowMap->getLocalNumElements(); ++i) {
490 size_t nnz =
tpOp.getNumEntriesInLocalRow(i);
491 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::nonconst_values_host_view_type
scaledVals(
"ScaledVals",nnz);
493 for (
size_t j = 0; j < nnz; ++j) {
503 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::global_inds_host_view_type
cols;
504 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::values_host_view_type
vals;
506 for (
size_t i = 0; i < rowMap->getLocalNumElements(); ++i) {
509 size_t nnz =
tpOp.getNumEntriesInGlobalRow(
gid);
510 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::nonconst_values_host_view_type
scaledVals(
"ScaledVals",nnz);
513 for (
size_t j = 0; j < nnz; ++j)
523 if (domainMap == Teuchos::null ||
rangeMap == Teuchos::null)
524 throw Exceptions::RuntimeError(
"In Utilities::Scaling: cannot fillComplete because the domain and/or range map hasn't been defined");
528 params->set(
"No Nonlocal Changes",
true);
529 Op.fillComplete(Op.getDomainMap(), Op.getRangeMap(),
params);
539 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
540 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
542 Transpose (Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op,
bool ,
const std::string &
label,
const Teuchos::RCP<Teuchos::ParameterList> &
params) {
543#if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT)
544 std::string
TorE =
"epetra";
546 std::string
TorE =
"tpetra";
549#if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT)
558#ifdef HAVE_MUELU_TPETRA
559 if (
TorE ==
"tpetra") {
560 using Helpers = Xpetra::Helpers<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
562 if(Helpers::isTpetraCrs(Op)) {
569 using Teuchos::ParameterList;
572 rcp (
new ParameterList) :
573 rcp (
new ParameterList (*
params));
581 if (!
AAAA->isFillComplete())
582 AAAA->fillComplete(Op.getRangeMap(), Op.getDomainMap());
584 if (Op.IsView(
"stridedMaps"))
585 AAAA->CreateView(
"stridedMaps", Teuchos::rcpFromRef(Op),
true);
589 else if(Helpers::isTpetraBlockCrs(Op)) {
590 using XMatrix = Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
591 using XCrsMatrix = Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
592 using XCrsMatrixWrap = Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
593 using BCRS = Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
594 using CRS = Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
599 Tpetra::BlockCrsMatrixTransposer<Scalar, LocalOrdinal, GlobalOrdinal, Node>
transposer(rcpFromRef(
tpetraOp),
label);
601 using Teuchos::ParameterList;
604 rcp (
new ParameterList) :
605 rcp (
new ParameterList (*
params));
614 if (Op.IsView(
"stridedMaps"))
615 AAAA->CreateView(
"stridedMaps", Teuchos::rcpFromRef(Op),
true);
625 std::cout <<
"Utilities::Transpose() not implemented for Epetra" << std::endl;
626 return Teuchos::null;
631 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
632 RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
636#if defined(HAVE_XPETRA_TPETRA) && (defined(HAVE_TPETRA_INST_COMPLEX_DOUBLE) || defined(HAVE_TPETRA_INST_COMPLEX_FLOAT))
637 typedef typename Teuchos::ScalarTraits<Scalar>::coordinateType real_type;
639 if ((
typeid(
Scalar).name() ==
typeid(std::complex<double>).name()) ||
640 (
typeid(
Scalar).name() ==
typeid(std::complex<float>).name())) {
641 Xscalar = Xpetra::MultiVectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(X->getMap(),X->getNumVectors());
642 size_t numVecs = X->getNumVectors();
643 for (
size_t j=0;j<
numVecs;j++) {
644 Teuchos::ArrayRCP<const real_type>
XVec = X->getData(j);
656 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
664 if(
paramList.isParameter (
"Coordinates") ==
false)
667#if defined(HAVE_MUELU_TPETRA)
673#if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT)
674 typedef Tpetra::MultiVector<float, LocalOrdinal, GlobalOrdinal, Node>
tfMV;
681#if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_DOUBLE)
689#if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT)
708 throw Exceptions::RuntimeError(
"ExtractCoordinatesFromParameterList: The coordinates vector in parameter list is expected to be a Tpetra multivector with SC=double or float.");
722#define MUELU_UTILITIES_SHORT
MueLu::DefaultLocalOrdinal LocalOrdinal
MueLu::DefaultScalar Scalar
MueLu::DefaultGlobalOrdinal GlobalOrdinal
Exception indicating invalid cast attempted.
Exception throws to report errors in the internal logical of the program.
Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > Map
static RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraCrs(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraBlockCrs(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraRow(RCP< Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static void MyOldScaleMatrix_Tpetra(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
static RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Transpose(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, bool optimizeTranspose=false, const std::string &label=std::string(), const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
static RCP< const Epetra_CrsMatrix > Op2EpetraCrs(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV2(Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
static const Epetra_Map & Map2EpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
static RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > RealValuedToScalarMultiVector(RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::coordinateType, LocalOrdinal, GlobalOrdinal, Node > > X)
static const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > Map2TpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
static RCP< const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2TpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > const vec)
Helper utility to pull out the underlying Tpetra objects from an Xpetra object.
static void MyOldScaleMatrix(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< const Scalar > &scalingVector, bool doInverse=true, bool doFillComplete=true, bool doOptimizeStorage=true)
static void MyOldScaleMatrix_Epetra(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec)
static RCP< Epetra_MultiVector > MV2NonConstEpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec)
static RCP< Epetra_CrsMatrix > Op2NonConstEpetraCrs(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraCrs(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraBlockCrs(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraRow(RCP< const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::magnitudeType, LocalOrdinal, GlobalOrdinal, Node > > ExtractCoordinatesFromParameterList(ParameterList ¶mList)
static RCP< const Epetra_MultiVector > MV2EpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > const vec)
Helper utility to pull out the underlying Epetra objects from an Xpetra object.
Namespace for MueLu classes and methods.
RCP< Xpetra::CrsMatrixWrap< SC, LO, GO, NO > > Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap(RCP< Epetra_CrsMatrix > &epAB)