46#ifndef MUELU_UTILITIES_DECL_HPP
47#define MUELU_UTILITIES_DECL_HPP
53#include <Teuchos_DefaultComm.hpp>
54#include <Teuchos_ScalarTraits.hpp>
55#include <Teuchos_ParameterList.hpp>
57#include <Xpetra_TpetraBlockCrsMatrix_fwd.hpp>
58#include <Xpetra_TpetraOperator.hpp>
59#include <Xpetra_CrsMatrix_fwd.hpp>
60#include <Xpetra_CrsMatrixWrap.hpp>
61#include <Xpetra_Map_fwd.hpp>
62#include <Xpetra_Matrix_fwd.hpp>
63#include <Xpetra_MultiVector_fwd.hpp>
64#include <Xpetra_MultiVectorFactory_fwd.hpp>
65#include <Xpetra_Operator_fwd.hpp>
66#include <Xpetra_Vector_fwd.hpp>
67#include <Xpetra_VectorFactory_fwd.hpp>
69#include <Xpetra_MatrixMatrix.hpp>
71#ifdef HAVE_MUELU_EPETRA
72#include <Xpetra_EpetraCrsMatrix.hpp>
76#include <Xpetra_EpetraCrsMatrix_fwd.hpp>
77#include <Xpetra_CrsMatrixWrap_fwd.hpp>
83#ifdef HAVE_MUELU_EPETRAEXT
90#include <Tpetra_CrsMatrix.hpp>
91#include <Tpetra_BlockCrsMatrix.hpp>
92#include <Tpetra_BlockCrsMatrix_Helpers.hpp>
93#include <Tpetra_FECrsMatrix.hpp>
94#include <Tpetra_RowMatrixTransposer.hpp>
95#include <Tpetra_Map.hpp>
96#include <Tpetra_MultiVector.hpp>
97#include <Tpetra_FEMultiVector.hpp>
98#include <Xpetra_TpetraCrsMatrix_fwd.hpp>
99#include <Xpetra_TpetraMultiVector_fwd.hpp>
101#include <MueLu_UtilitiesBase.hpp>
106#ifdef HAVE_MUELU_EPETRA
108 template<
typename SC,
typename LO,
typename GO,
typename NO>
109 RCP<Xpetra::CrsMatrixWrap<SC,LO,GO,NO> >
112 template<
typename SC,
typename LO,
typename GO,
typename NO>
113 RCP<Xpetra::Matrix<SC, LO, GO, NO> >
116 template<
typename SC,
typename LO,
typename GO,
typename NO>
117 RCP<Xpetra::MultiVector<SC, LO, GO, NO> >
121 template<
typename SC,
typename LO,
typename GO,
typename NO>
122 RCP<Xpetra::Matrix<SC, LO, GO, NO> >
125 template<
typename SC,
typename LO,
typename GO,
typename NO>
126 RCP<Xpetra::Matrix<SC, LO, GO, NO> >
129 template<
typename SC,
typename LO,
typename GO,
typename NO>
130 RCP<Xpetra::MultiVector<SC, LO, GO, NO> >
133 template<
typename SC,
typename LO,
typename GO,
typename NO>
134 RCP<Xpetra::MultiVector<SC, LO, GO, NO> >
137 template<
typename SC,
typename LO,
typename GO,
typename NO>
138 void leftRghtDofScalingWithinNode(
const Xpetra::Matrix<SC,LO,GO,NO> & Atpetra,
size_t blkSize,
size_t nSweeps, Teuchos::ArrayRCP<SC> & rowScaling, Teuchos::ArrayRCP<SC> & colScaling);
152#undef MUELU_UTILITIES_SHORT
156 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType
Magnitude;
158#ifdef HAVE_MUELU_EPETRA
161 static RCP<const Epetra_MultiVector>
MV2EpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
const vec);
162 static RCP< Epetra_MultiVector>
MV2NonConstEpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > vec);
167 static RCP<const Epetra_CrsMatrix>
Op2EpetraCrs(RCP<
const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
168 static RCP< Epetra_CrsMatrix>
Op2NonConstEpetraCrs(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
178 static RCP<const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
MV2TpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
const vec);
179 static RCP< Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
MV2NonConstTpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > vec);
180 static RCP< Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
MV2NonConstTpetraMV2(Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
182 static const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
MV2TpetraMV(
const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
183 static Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
MV2NonConstTpetraMV(Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
185 static RCP<const Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
Op2TpetraCrs(RCP<
const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
186 static RCP< Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
Op2NonConstTpetraCrs(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
188 static const Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
Op2TpetraCrs(
const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
189 static Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
Op2NonConstTpetraCrs(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
191 static RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
Op2TpetraBlockCrs(RCP<
const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
192 static RCP< Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
Op2NonConstTpetraBlockCrs(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
194 static const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
Op2TpetraBlockCrs(
const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
195 static Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
Op2NonConstTpetraBlockCrs(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
199 static RCP<const Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
Op2TpetraRow(RCP<
const Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
200 static RCP< Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
Op2NonConstTpetraRow(RCP<Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
203 static const RCP<const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> >
Map2TpetraMap(
const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node>& map);
205 static void MyOldScaleMatrix(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op,
const Teuchos::ArrayRCP<const Scalar>& scalingVector,
bool doInverse =
true,
206 bool doFillComplete =
true,
bool doOptimizeStorage =
true);
208 static void MyOldScaleMatrix_Epetra(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op,
const Teuchos::ArrayRCP<Scalar>& scalingVector,
209 bool doFillComplete,
bool doOptimizeStorage);
210 static void MyOldScaleMatrix_Tpetra(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op,
const Teuchos::ArrayRCP<Scalar>& scalingVector,
211 bool doFillComplete,
bool doOptimizeStorage);
213 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);
223#ifdef HAVE_MUELU_EPETRA
239 typedef Xpetra::EpetraNode
Node;
240 typedef Teuchos::ScalarTraits<Scalar>::magnitudeType
Magnitude;
241#undef MUELU_UTILITIES_SHORT
245 using EpetraMap = Xpetra::EpetraMapT<GlobalOrdinal,Node>;
252 static RCP<const Epetra_MultiVector>
MV2EpetraMV(RCP<MultiVector>
const vec) {
253 RCP<const EpetraMultiVector > tmpVec = rcp_dynamic_cast<EpetraMultiVector>(vec);
254 if (tmpVec == Teuchos::null)
255 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
256 return tmpVec->getEpetra_MultiVector();
259 RCP<const EpetraMultiVector> tmpVec = rcp_dynamic_cast<EpetraMultiVector>(vec);
260 if (tmpVec == Teuchos::null)
261 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
262 return tmpVec->getEpetra_MultiVector();
267 return *(tmpVec.getEpetra_MultiVector());
271 return *(tmpVec.getEpetra_MultiVector());
274 static RCP<const Epetra_CrsMatrix>
Op2EpetraCrs(RCP<const Matrix> Op) {
275 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
276 if (crsOp == Teuchos::null)
278 const RCP<const EpetraCrsMatrix>& tmp_ECrsMtx = rcp_dynamic_cast<const EpetraCrsMatrix>(crsOp->getCrsMatrix());
279 if (tmp_ECrsMtx == Teuchos::null)
280 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
281 return tmp_ECrsMtx->getEpetra_CrsMatrix();
284 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
285 if (crsOp == Teuchos::null)
287 const RCP<const EpetraCrsMatrix> &tmp_ECrsMtx = rcp_dynamic_cast<const EpetraCrsMatrix>(crsOp->getCrsMatrix());
288 if (tmp_ECrsMtx == Teuchos::null)
289 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
290 return tmp_ECrsMtx->getEpetra_CrsMatrixNonConst();
295 const CrsMatrixWrap& crsOp =
dynamic_cast<const CrsMatrixWrap&
>(Op);
297 const EpetraCrsMatrix& tmp_ECrsMtx =
dynamic_cast<const EpetraCrsMatrix&
>(*crsOp.getCrsMatrix());
298 return *tmp_ECrsMtx.getEpetra_CrsMatrix();
299 }
catch (std::bad_cast&) {
300 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
302 }
catch (std::bad_cast&) {
308 CrsMatrixWrap& crsOp =
dynamic_cast<CrsMatrixWrap&
>(Op);
310 EpetraCrsMatrix& tmp_ECrsMtx =
dynamic_cast<EpetraCrsMatrix&
>(*crsOp.getCrsMatrix());
311 return *tmp_ECrsMtx.getEpetra_CrsMatrixNonConst();
312 }
catch (std::bad_cast&) {
313 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
315 }
catch (std::bad_cast&) {
321 RCP<const EpetraMap> xeMap = rcp_dynamic_cast<const EpetraMap>(rcpFromRef(map));
322 if (xeMap == Teuchos::null)
323 throw Exceptions::BadCast(
"Utilities::Map2EpetraMap : Cast from Xpetra::Map to Xpetra::EpetraMap failed");
324 return xeMap->getEpetra_Map();
329 static RCP<const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
MV2TpetraMV(RCP<MultiVector>
const vec) {
330#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
331 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
334 RCP<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(vec);
335 if (tmpVec == Teuchos::null)
336 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
337 return tmpVec->getTpetra_MultiVector();
340 static RCP< Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
MV2NonConstTpetraMV(RCP<MultiVector> vec) {
341#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
342 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
345 RCP<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(vec);
346 if (tmpVec == Teuchos::null)
347 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
348 return tmpVec->getTpetra_MultiVector();
352 static RCP< Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
MV2NonConstTpetraMV2(MultiVector& vec) {
353#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
354 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
357 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec =
dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(vec);
358 return tmpVec.getTpetra_MultiVector();
362 static const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
MV2TpetraMV(
const MultiVector& vec) {
363#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
364 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
367 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec =
dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(vec);
368 return *(tmpVec.getTpetra_MultiVector());
371 static Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
MV2NonConstTpetraMV(MultiVector& vec) {
372#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
373 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
376 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec =
dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(vec);
377 return *(tmpVec.getTpetra_MultiVector());
381 static RCP<const Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
Op2TpetraCrs(RCP<const Matrix> Op) {
382#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
383 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
387 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
388 if (crsOp == Teuchos::null)
390 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
391 if (tmp_ECrsMtx == Teuchos::null)
392 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
393 return tmp_ECrsMtx->getTpetra_CrsMatrix();
397#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
398 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
401 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
402 if (crsOp == Teuchos::null)
404 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
405 if (tmp_ECrsMtx == Teuchos::null)
406 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
407 return tmp_ECrsMtx->getTpetra_CrsMatrixNonConst();
411 static const Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
Op2TpetraCrs(
const Matrix& Op) {
412#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
413 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
417 const CrsMatrixWrap& crsOp =
dynamic_cast<const CrsMatrixWrap&
>(Op);
419 const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx =
dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(*crsOp.getCrsMatrix());
420 return *tmp_ECrsMtx.getTpetra_CrsMatrix();
421 }
catch (std::bad_cast&) {
422 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
424 }
catch (std::bad_cast&) {
430#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
431 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
435 CrsMatrixWrap& crsOp =
dynamic_cast<CrsMatrixWrap&
>(Op);
437 Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx =
dynamic_cast<Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(*crsOp.getCrsMatrix());
438 return *tmp_ECrsMtx.getTpetra_CrsMatrixNonConst();
439 }
catch (std::bad_cast&) {
440 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
442 }
catch (std::bad_cast&) {
449 static RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
Op2TpetraBlockCrs(RCP<const Matrix> Op) {
450#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
451 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
455 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
456 if (crsOp == Teuchos::null)
458 const RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
459 if (tmp_ECrsMtx == Teuchos::null)
460 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
461 return tmp_ECrsMtx->getTpetra_BlockCrsMatrix();
466#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
467 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
468 throw Exceptions::RuntimeError(
"Op2NonConstTpetraBlockCrs: Tpetra has not been compiled with support for LO=GO=int.");
470 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
471 if (crsOp == Teuchos::null)
473 const RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
474 if (tmp_ECrsMtx == Teuchos::null)
475 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
476 return tmp_ECrsMtx->getTpetra_BlockCrsMatrixNonConst();
480 static const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
Op2TpetraBlockCrs(
const Matrix& Op) {
481#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
482 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
486 const CrsMatrixWrap& crsOp =
dynamic_cast<const CrsMatrixWrap&
>(Op);
488 const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx =
dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(*crsOp.getCrsMatrix());
489 return *tmp_ECrsMtx.getTpetra_BlockCrsMatrix();
490 }
catch (std::bad_cast&) {
491 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
493 }
catch (std::bad_cast&) {
499#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
500 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
504 CrsMatrixWrap& crsOp =
dynamic_cast<CrsMatrixWrap&
>(Op);
506 Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx =
dynamic_cast<Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(*crsOp.getCrsMatrix());
507 return *tmp_ECrsMtx.getTpetra_BlockCrsMatrixNonConst();
508 }
catch (std::bad_cast&) {
509 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
511 }
catch (std::bad_cast&) {
518 static RCP<const Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
Op2TpetraRow(RCP<const Operator> Op) {
519#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
520 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
523 RCP<const Matrix> mat = rcp_dynamic_cast<const Matrix>(Op);
524 if (!mat.is_null()) {
525 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(mat);
526 if (crsOp == Teuchos::null)
529 RCP<const CrsMatrix> crsMat = crsOp->getCrsMatrix();
530 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_Crs = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
531 RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_BlockCrs;
532 if(!tmp_Crs.is_null()) {
533 return tmp_Crs->getTpetra_CrsMatrixNonConst();
536 tmp_BlockCrs= rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
537 if (tmp_BlockCrs.is_null())
538 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
539 return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
542 RCP<const TpetraOperator> tpOp = rcp_dynamic_cast<const TpetraOperator>(Op,
true);
543 RCP<const Tpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node>> tOp = tpOp->getOperatorConst();
544 RCP<const Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tRow = rcp_dynamic_cast<const Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(tOp,
true);
550 static RCP< Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
Op2NonConstTpetraRow(RCP<Operator> Op) {
551#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
552 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
555 RCP<Matrix> mat = rcp_dynamic_cast<Matrix>(Op);
556 if (!mat.is_null()) {
557 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(mat);
558 if (crsOp == Teuchos::null)
561 RCP<const CrsMatrix> crsMat = crsOp->getCrsMatrix();
562 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_Crs = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
563 RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_BlockCrs;
564 if(!tmp_Crs.is_null()) {
565 return tmp_Crs->getTpetra_CrsMatrixNonConst();
568 tmp_BlockCrs= rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
569 if (tmp_BlockCrs.is_null())
570 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
571 return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
574 RCP<TpetraOperator> tpOp = rcp_dynamic_cast<TpetraOperator>(Op,
true);
575 RCP<Tpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node>> tOp = tpOp->getOperator();
576 RCP<Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tRow = rcp_dynamic_cast<Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(tOp,
true);
583 static const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> >
Map2TpetraMap(
const Map& map) {
584#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
585 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
588 const RCP<const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node>>& tmp_TMap = rcp_dynamic_cast<const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node> >(rcpFromRef(map));
589 if (tmp_TMap == Teuchos::null)
590 throw Exceptions::BadCast(
"Utilities::Map2TpetraMap : Cast from Xpetra::Map to Xpetra::TpetraMap failed");
591 return tmp_TMap->getTpetra_Map();
596 static void MyOldScaleMatrix(Matrix& Op,
const Teuchos::ArrayRCP<const Scalar>& scalingVector,
bool doInverse =
true,
597 bool doFillComplete =
true,
bool doOptimizeStorage =
true) {
598 Scalar one = Teuchos::ScalarTraits<Scalar>::one();
599 Teuchos::ArrayRCP<Scalar> sv(scalingVector.size());
601 for (
int i = 0; i < scalingVector.size(); ++i)
602 sv[i] = one / scalingVector[i];
604 for (
int i = 0; i < scalingVector.size(); ++i)
605 sv[i] = scalingVector[i];
608 switch (Op.getRowMap()->lib()) {
609 case Xpetra::UseTpetra:
613 case Xpetra::UseEpetra:
624 bool doFillComplete,
bool doOptimizeStorage) {
625#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
626 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
627 throw Exceptions::RuntimeError(
"Matrix scaling is not possible because Tpetra has not been compiled with support for LO=GO=int.");
632 const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rowMap = tpOp.getRowMap();
633 const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > domainMap = tpOp.getDomainMap();
634 const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rangeMap = tpOp.getRangeMap();
636 size_t maxRowSize = tpOp.getLocalMaxNumRowEntries();
637 if (maxRowSize == Teuchos::as<size_t>(-1))
640 std::vector<Scalar> scaledVals(maxRowSize);
641 if (tpOp.isFillComplete())
644 if (Op.isLocallyIndexed() ==
true) {
645 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::local_inds_host_view_type cols;
646 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::values_host_view_type vals;
647 for (
size_t i = 0; i < rowMap->getLocalNumElements(); ++i) {
648 tpOp.getLocalRowView(i, cols, vals);
649 size_t nnz = tpOp.getNumEntriesInLocalRow(i);
650 if (nnz > maxRowSize) {
652 scaledVals.resize(maxRowSize);
654 for (
size_t j = 0; j < nnz; ++j)
655 scaledVals[j] = vals[j]*scalingVector[i];
658 Teuchos::ArrayView<const LocalOrdinal> cols_view(cols.data(), nnz);
659 Teuchos::ArrayView<const Scalar> valview(&scaledVals[0], nnz);
660 tpOp.replaceLocalValues(i, cols_view, valview);
665 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::local_inds_host_view_type cols;
666 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::values_host_view_type vals;
668 for (
size_t i = 0; i < rowMap->getLocalNumElements(); ++i) {
670 tpOp.getGlobalRowView(gid, cols, vals);
671 size_t nnz = tpOp.getNumEntriesInGlobalRow(gid);
672 if (nnz > maxRowSize) {
674 scaledVals.resize(maxRowSize);
677 for (
size_t j = 0; j < nnz; ++j)
678 scaledVals[j] = vals[j]*scalingVector[i];
681 Teuchos::ArrayView<const LocalOrdinal> cols_view(cols.data(), nnz);
682 Teuchos::ArrayView<const Scalar> valview(&scaledVals[0], nnz);
683 tpOp.replaceGlobalValues(gid, cols_view, valview);
688 if (doFillComplete) {
689 if (domainMap == Teuchos::null || rangeMap == Teuchos::null)
690 throw Exceptions::RuntimeError(
"In Utilities::Scaling: cannot fillComplete because the domain and/or range map hasn't been defined");
692 RCP<Teuchos::ParameterList> params = rcp(
new Teuchos::ParameterList());
693 params->set(
"Optimize Storage", doOptimizeStorage);
694 params->set(
"No Nonlocal Changes",
true);
695 Op.fillComplete(Op.getDomainMap(), Op.getRangeMap(), params);
704#ifdef HAVE_MUELU_EPETRA
716 for (
int j = 0; j < nnz; ++j)
717 vals[j] *= scalingVector[i];
733 static RCP<Matrix>
Transpose(Matrix& Op,
bool =
false,
const std::string & label = std::string(),
const Teuchos::RCP<Teuchos::ParameterList> ¶ms=Teuchos::null) {
734 switch (Op.getRowMap()->lib()) {
735 case Xpetra::UseTpetra: {
736#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
737 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
738 throw Exceptions::RuntimeError(
"Utilities::Transpose: Tpetra is not compiled with LO=GO=int. Add TPETRA_INST_INT_INT:BOOL=ON to your configuration!");
740 using Helpers = Xpetra::Helpers<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
742 if(Helpers::isTpetraCrs(Op)) {
746 RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > A;
747 Tpetra::RowMatrixTransposer<Scalar, LocalOrdinal, GlobalOrdinal, Node> transposer(rcpFromRef(tpetraOp),label);
750 using Teuchos::ParameterList;
752 RCP<ParameterList> transposeParams = params.is_null () ?
753 rcp (
new ParameterList) :
754 rcp (
new ParameterList (*params));
755 transposeParams->set (
"sort",
false);
756 A = transposer.createTranspose(transposeParams);
759 RCP<Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AA = rcp(
new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(A));
760 RCP<CrsMatrix> AAA = rcp_implicit_cast<CrsMatrix>(AA);
761 RCP<Matrix> AAAA = rcp(
new CrsMatrixWrap(AAA));
763 if (Op.IsView(
"stridedMaps"))
764 AAAA->CreateView(
"stridedMaps", Teuchos::rcpFromRef(Op),
true);
769 else if(Helpers::isTpetraBlockCrs(Op)) {
770 using BCRS = Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
775 Tpetra::BlockCrsMatrixTransposer<Scalar, LocalOrdinal, GlobalOrdinal, Node> transposer(rcpFromRef(tpetraOp),label);
777 using Teuchos::ParameterList;
779 RCP<ParameterList> transposeParams = params.is_null () ?
780 rcp (
new ParameterList) :
781 rcp (
new ParameterList (*params));
782 transposeParams->set (
"sort",
false);
783 At = transposer.createTranspose(transposeParams);
786 RCP<Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AA = rcp(
new Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(At));
787 RCP<CrsMatrix> AAA = rcp_implicit_cast<CrsMatrix>(AA);
788 RCP<Matrix> AAAA = rcp(
new CrsMatrixWrap(AAA));
790 if (Op.IsView(
"stridedMaps"))
791 AAAA->CreateView(
"stridedMaps", Teuchos::rcpFromRef(Op),
true);
798 throw Exceptions::RuntimeError(
"Utilities::Transpose failed, perhaps because matrix is not a Crs or BlockCrs matrix");
802 case Xpetra::UseEpetra:
804#if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT)
805 Teuchos::TimeMonitor tm(*Teuchos::TimeMonitor::getNewTimer(
"ZZ Entire Transpose"));
812 RCP<Epetra_CrsMatrix> rcpA(A);
813 RCP<EpetraCrsMatrix> AA = rcp(
new EpetraCrsMatrix(rcpA));
814 RCP<CrsMatrix> AAA = rcp_implicit_cast<CrsMatrix>(AA);
815 RCP<Matrix> AAAA = rcp(
new CrsMatrixWrap(AAA));
817 if (Op.IsView(
"stridedMaps"))
818 AAAA->CreateView(
"stridedMaps", Teuchos::rcpFromRef(Op),
true);
829 TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
834 RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Xscalar = rcp_dynamic_cast<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(X,
true);
841 RCP<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,
LocalOrdinal,
GlobalOrdinal,
Node> > coordinates = Teuchos::null;
844 if(paramList.isParameter (
"Coordinates") ==
false)
847 #if ( defined(EPETRA_HAVE_OMP) && defined(HAVE_TPETRA_INST_OPENMP) && defined(HAVE_TPETRA_INST_INT_INT)) || \
848 (!defined(EPETRA_HAVE_OMP) && defined(HAVE_TPETRA_INST_SERIAL) && defined(HAVE_TPETRA_INST_INT_INT))
853 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT)
854 typedef Tpetra::MultiVector<float, LocalOrdinal, GlobalOrdinal, Node> tfMV;
855 RCP<tfMV> floatCoords = Teuchos::null;
862 RCP<tdMV> doubleCoords = Teuchos::null;
863 if (paramList.isType<RCP<tdMV> >(
"Coordinates")) {
865 doubleCoords = paramList.get<RCP<tdMV> >(
"Coordinates");
866 paramList.remove(
"Coordinates");
868 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT)
869 else if (paramList.isType<RCP<tfMV> >(
"Coordinates")) {
871 floatCoords = paramList.get<RCP<tfMV> >(
"Coordinates");
872 paramList.remove(
"Coordinates");
873 doubleCoords = rcp(
new tdMV(floatCoords->getMap(), floatCoords->getNumVectors()));
874 deep_copy(*doubleCoords, *floatCoords);
878 if(doubleCoords != Teuchos::null) {
879 coordinates = Teuchos::rcp(
new Xpetra::TpetraMultiVector<
typename Teuchos::ScalarTraits<Scalar>::magnitudeType,
LocalOrdinal,
GlobalOrdinal,
Node>(doubleCoords));
880 TEUCHOS_TEST_FOR_EXCEPT(doubleCoords->getNumVectors() != coordinates->getNumVectors());
884 #if defined(HAVE_MUELU_EPETRA)
885 RCP<Epetra_MultiVector> doubleEpCoords;
886 if (paramList.isType<RCP<Epetra_MultiVector> >(
"Coordinates")) {
887 doubleEpCoords = paramList.get<RCP<Epetra_MultiVector> >(
"Coordinates");
888 paramList.remove(
"Coordinates");
889 RCP<Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> > epCoordinates = Teuchos::rcp(
new Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node>(doubleEpCoords));
890 coordinates = rcp_dynamic_cast<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,
LocalOrdinal,
GlobalOrdinal,
Node> >(epCoordinates);
891 TEUCHOS_TEST_FOR_EXCEPT(doubleEpCoords->NumVectors() != Teuchos::as<int>(coordinates->getNumVectors()));
896 if(paramList.isType<
decltype(coordinates)>(
"Coordinates")) {
897 coordinates = paramList.get<
decltype(coordinates)>(
"Coordinates");
900 TEUCHOS_TEST_FOR_EXCEPT(Teuchos::is_null(coordinates));
939 long ExtractNonSerializableData(
const Teuchos::ParameterList& inList, Teuchos::ParameterList& serialList, Teuchos::ParameterList& nonSerialList);
955#ifdef HAVE_MUELU_EPETRA
960 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
961 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
963 typedef Xpetra::EpetraCrsMatrixT<GlobalOrdinal, Node> XECrsMatrix;
964 typedef Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrix;
965 typedef Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrixWrap;
967 RCP<XCrsMatrix> Atmp = rcp(
new XECrsMatrix(A));
968 return rcp(
new XCrsMatrixWrap(Atmp));
975 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
976 RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
978 return rcp(
new Xpetra::EpetraMultiVectorT<GlobalOrdinal, Node>(V));
986 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
987 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
989 typedef Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> XTCrsMatrix;
990 typedef Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrix;
991 typedef Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrixWrap;
993 RCP<XCrsMatrix> Atmp = rcp(
new XTCrsMatrix(Atpetra));
994 return rcp(
new XCrsMatrixWrap(Atmp));
1024 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1025 void leftRghtDofScalingWithinNode(
const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> & Amat,
size_t blkSize,
size_t nSweeps, Teuchos::ArrayRCP<Scalar> & rowScaling, Teuchos::ArrayRCP<Scalar> & colScaling) {
1027 LocalOrdinal nBlks = (Amat.getRowMap()->getLocalNumElements())/blkSize;
1029 Teuchos::ArrayRCP<Scalar> rowScaleUpdate(blkSize);
1030 Teuchos::ArrayRCP<Scalar> colScaleUpdate(blkSize);
1033 for (
size_t i = 0; i < blkSize; i++) rowScaling[i] = 1.0;
1034 for (
size_t i = 0; i < blkSize; i++) colScaling[i] = 1.0;
1036 for (
size_t k = 0; k < nSweeps; k++) {
1038 for (
size_t i = 0; i < blkSize; i++) rowScaleUpdate[i] = 0.0;
1041 for (
size_t j = 0; j < blkSize; j++) {
1042 Teuchos::ArrayView<const LocalOrdinal> cols;
1043 Teuchos::ArrayView<const Scalar> vals;
1044 Amat.getLocalRowView(row, cols, vals);
1046 for (
size_t kk = 0; kk < Teuchos::as<size_t>(vals.size()); kk++) {
1047 size_t modGuy = (cols[kk]+1)%blkSize;
1048 if (modGuy == 0) modGuy = blkSize;
1050 rowScaleUpdate[j] += rowScaling[j]*(Teuchos::ScalarTraits<Scalar>::magnitude(vals[kk]))*colScaling[modGuy];
1056 Teuchos::ArrayRCP<Scalar> tempUpdate(blkSize);
1057 Teuchos::reduceAll(*(Amat.getRowMap()->getComm()), Teuchos::REDUCE_SUM, (
LocalOrdinal) blkSize, rowScaleUpdate.getRawPtr(), tempUpdate.getRawPtr());
1058 for (
size_t i = 0; i < blkSize; i++) rowScaleUpdate[i] = tempUpdate[i];
1064 Scalar minUpdate = Teuchos::ScalarTraits<Scalar>::magnitude((rowScaleUpdate[0]/rowScaling[0])/rowScaling[0]);
1066 for (
size_t i = 1; i < blkSize; i++) {
1067 Scalar temp = (rowScaleUpdate[i]/rowScaling[i])/rowScaling[i];
1068 if ( Teuchos::ScalarTraits<Scalar>::magnitude(temp) < Teuchos::ScalarTraits<Scalar>::magnitude(minUpdate))
1069 minUpdate = Teuchos::ScalarTraits<Scalar>::magnitude(temp);
1071 for (
size_t i = 0; i < blkSize; i++) rowScaling[i] *= sqrt(minUpdate / rowScaleUpdate[i]);
1074 for (
size_t i = 0; i < blkSize; i++) colScaleUpdate[i] = 0.0;
1077 for (
size_t j = 0; j < blkSize; j++) {
1078 Teuchos::ArrayView<const LocalOrdinal> cols;
1079 Teuchos::ArrayView<const Scalar> vals;
1080 Amat.getLocalRowView(row, cols, vals);
1081 for (
size_t kk = 0; kk < Teuchos::as<size_t>(vals.size()); kk++) {
1082 size_t modGuy = (cols[kk]+1)%blkSize;
1083 if (modGuy == 0) modGuy = blkSize;
1085 colScaleUpdate[modGuy] += colScaling[modGuy]* (Teuchos::ScalarTraits<Scalar>::magnitude(vals[kk])) *rowScaling[j];
1090 Teuchos::reduceAll(*(Amat.getRowMap()->getComm()), Teuchos::REDUCE_SUM, (
LocalOrdinal) blkSize, colScaleUpdate.getRawPtr(), tempUpdate.getRawPtr());
1091 for (
size_t i = 0; i < blkSize; i++) colScaleUpdate[i] = tempUpdate[i];
1098 minUpdate = Teuchos::ScalarTraits<Scalar>::magnitude((colScaleUpdate[0]/colScaling[0])/colScaling[0]);
1100 for (
size_t i = 1; i < blkSize; i++) {
1101 Scalar temp = (colScaleUpdate[i]/colScaling[i])/colScaling[i];
1102 if ( Teuchos::ScalarTraits<Scalar>::magnitude(temp) < Teuchos::ScalarTraits<Scalar>::magnitude(minUpdate))
1103 minUpdate = Teuchos::ScalarTraits<Scalar>::magnitude(temp);
1105 for (
size_t i = 0; i < blkSize; i++) colScaling[i] *= sqrt(minUpdate/colScaleUpdate[i]);
1113 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1114 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1116 typedef typename Tpetra::FECrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::crs_matrix_type tpetra_crs_matrix_type;
1117 typedef Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> XTCrsMatrix;
1118 typedef Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrix;
1119 typedef Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrixWrap;
1121 RCP<XCrsMatrix> Atmp = rcp(
new XTCrsMatrix(rcp_dynamic_cast<tpetra_crs_matrix_type>(Atpetra)));
1122 return rcp(
new XCrsMatrixWrap(Atmp));
1129 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1130 RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1132 return rcp(
new Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>(Vtpetra));
1139 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1140 RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1142 typedef Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> MV;
1143 RCP<const MV> Vmv = Teuchos::rcp_dynamic_cast<const MV>(Vtpetra);
1144 return rcp(
new Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>(Vmv));
1150 std::ostringstream buf;
1155#ifdef HAVE_MUELU_EPETRA
1160 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1161 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1168 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1169 RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1177 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1178 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1185 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1186 RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1193 std::string
lowerCase (
const std::string& s);
1197#define MUELU_UTILITIES_SHORT
Tpetra::KokkosCompat::KokkosSerialWrapperNode EpetraNode
MueLu::DefaultLocalOrdinal LocalOrdinal
MueLu::DefaultScalar Scalar
MueLu::DefaultGlobalOrdinal GlobalOrdinal
int NumMyElements() const
const Epetra_Map & RowMap() const
int ExtractMyRowView(int MyRow, int &NumEntries, double *&Values, int *&Indices) const
Exception indicating invalid cast attempted.
Exception throws to report errors in the internal logical of the program.
static const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & MV2TpetraMV(const MultiVector &vec)
static RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraBlockCrs(RCP< const Matrix > Op)
static RCP< Epetra_MultiVector > MV2NonConstEpetraMV(RCP< MultiVector > vec)
static RCP< Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraRow(RCP< Operator > Op)
static Epetra_CrsMatrix & Op2NonConstEpetraCrs(Matrix &Op)
Teuchos::ScalarTraits< Scalar >::magnitudeType Magnitude
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV2(MultiVector &vec)
static void MyOldScaleMatrix(Matrix &Op, const Teuchos::ArrayRCP< const Scalar > &scalingVector, bool doInverse=true, bool doFillComplete=true, bool doOptimizeStorage=true)
static Epetra_MultiVector & MV2NonConstEpetraMV(MultiVector &vec)
static void MyOldScaleMatrix_Epetra(Matrix &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool, bool)
static RCP< Matrix > Transpose(Matrix &Op, bool=false, const std::string &label=std::string(), const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Transpose a Xpetra::Matrix.
static void MyOldScaleMatrix_Tpetra(Matrix &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV(RCP< MultiVector > vec)
static RCP< const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraRow(RCP< const Operator > Op)
Xpetra::EpetraMapT< GlobalOrdinal, Node > EpetraMap
static Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2NonConstTpetraBlockCrs(Matrix &Op)
static RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraCrs(RCP< Matrix > Op)
static const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2TpetraBlockCrs(const Matrix &Op)
static RCP< const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2TpetraMV(RCP< MultiVector > const vec)
Helper utility to pull out the underlying Tpetra objects from an Xpetra object.
static const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2TpetraCrs(const Matrix &Op)
static RCP< const Epetra_MultiVector > MV2EpetraMV(RCP< MultiVector > const vec)
Helper utility to pull out the underlying Epetra objects from an Xpetra object.
static const Epetra_Map & Map2EpetraMap(const Map &map)
static const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > Map2TpetraMap(const Map &map)
static RCP< const Epetra_CrsMatrix > Op2EpetraCrs(RCP< const Matrix > Op)
static const Epetra_MultiVector & MV2EpetraMV(const MultiVector &vec)
static Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2NonConstTpetraCrs(Matrix &Op)
static RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraCrs(RCP< const Matrix > Op)
Xpetra::EpetraMultiVectorT< GlobalOrdinal, Node > EpetraMultiVector
static RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::magnitudeType, LocalOrdinal, GlobalOrdinal, Node > > ExtractCoordinatesFromParameterList(ParameterList ¶mList)
Extract coordinates from parameter list and return them in a Xpetra::MultiVector.
static RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::magnitudeType, LocalOrdinal, GlobalOrdinal, Node > > RealValuedToScalarMultiVector(RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::coordinateType, LocalOrdinal, GlobalOrdinal, Node > > X)
static RCP< Epetra_CrsMatrix > Op2NonConstEpetraCrs(RCP< Matrix > Op)
static Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & MV2NonConstTpetraMV(MultiVector &vec)
static const Epetra_CrsMatrix & Op2EpetraCrs(const Matrix &Op)
static RCP< Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraBlockCrs(RCP< Matrix > Op)
static const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2TpetraCrs(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op)
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)
Teuchos::ScalarTraits< Scalar >::magnitudeType Magnitude
static const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & MV2TpetraMV(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
static void MyOldScaleMatrix_Tpetra(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
static Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2NonConstTpetraCrs(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op)
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 const Epetra_MultiVector & MV2EpetraMV(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
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 Epetra_CrsMatrix & Op2NonConstEpetraCrs(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op)
static const Epetra_Map & Map2EpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
static Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2NonConstTpetraBlockCrs(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op)
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 Epetra_MultiVector & MV2NonConstEpetraMV(Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
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 Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & MV2NonConstTpetraMV(Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
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 const Epetra_CrsMatrix & Op2EpetraCrs(const Xpetra::Matrix< 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.
static const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2TpetraBlockCrs(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op)
Namespace for MueLu classes and methods.
bool IsParamMuemexVariable(const std::string &name)
long ExtractNonSerializableData(const Teuchos::ParameterList &inList, Teuchos::ParameterList &serialList, Teuchos::ParameterList &nonSerialList)
Extract non-serializable data from level-specific sublists and move it to a separate parameter list.
std::string lowerCase(const std::string &s)
RCP< Xpetra::MultiVector< SC, LO, GO, NO > > TpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP< Tpetra::MultiVector< SC, LO, GO, NO > > &Vtpetra)
void leftRghtDofScalingWithinNode(const Xpetra::Matrix< SC, LO, GO, NO > &Atpetra, size_t blkSize, size_t nSweeps, Teuchos::ArrayRCP< SC > &rowScaling, Teuchos::ArrayRCP< SC > &colScaling)
RCP< Xpetra::MultiVector< SC, LO, GO, NO > > TpetraFEMultiVector_To_XpetraMultiVector(const Teuchos::RCP< Tpetra::FEMultiVector< SC, LO, GO, NO > > &Vtpetra)
RCP< Xpetra::MultiVector< SC, LO, GO, NO > > EpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP< Epetra_MultiVector > &V)
Tpetra::KokkosClassic::DefaultNode::DefaultNodeType DefaultNode
RCP< Xpetra::Matrix< SC, LO, GO, NO > > EpetraCrs_To_XpetraMatrix(const Teuchos::RCP< Epetra_CrsMatrix > &A)
bool IsParamValidVariable(const std::string &name)
void TokenizeStringAndStripWhiteSpace(const std::string &stream, std::vector< std::string > &tokenList, const char *delimChars)
RCP< Xpetra::Matrix< SC, LO, GO, NO > > TpetraCrs_To_XpetraMatrix(const Teuchos::RCP< Tpetra::CrsMatrix< SC, LO, GO, NO > > &Atpetra)
RCP< Xpetra::Matrix< SC, LO, GO, NO > > TpetraFECrs_To_XpetraMatrix(const Teuchos::RCP< Tpetra::FECrsMatrix< SC, LO, GO, NO > > &Atpetra)
RCP< Xpetra::CrsMatrixWrap< SC, LO, GO, NO > > Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap(RCP< Epetra_CrsMatrix > &epAB)
std::string toString(const T &what)
Little helper function to convert non-string types to strings.
Teuchos::RCP< const Teuchos::Comm< int > > GenerateNodeComm(RCP< const Teuchos::Comm< int > > &baseComm, int &NodeId, const int reductionFactor)