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_TpetraCrsMatrix.hpp> 78 #include <Xpetra_TpetraBlockCrsMatrix.hpp> 81 #ifdef HAVE_MUELU_EPETRA 82 #include <Xpetra_EpetraMap.hpp> 85 #include <Xpetra_BlockedCrsMatrix.hpp> 87 #include <Xpetra_IO.hpp> 88 #include <Xpetra_Import.hpp> 89 #include <Xpetra_ImportFactory.hpp> 90 #include <Xpetra_Map.hpp> 91 #include <Xpetra_MapFactory.hpp> 92 #include <Xpetra_Matrix.hpp> 93 #include <Xpetra_MatrixFactory.hpp> 94 #include <Xpetra_MultiVector.hpp> 95 #include <Xpetra_MultiVectorFactory.hpp> 96 #include <Xpetra_Operator.hpp> 97 #include <Xpetra_Vector.hpp> 98 #include <Xpetra_VectorFactory.hpp> 100 #include <Xpetra_MatrixMatrix.hpp> 103 #if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_ML) 104 #include <ml_operator.h> 105 #include <ml_epetra_utils.h> 110 #ifdef HAVE_MUELU_EPETRA 115 #ifdef HAVE_MUELU_EPETRA 116 template<
typename SC,
typename LO,
typename GO,
typename NO>
118 return Xpetra::Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap<SC,LO,GO,NO>(epAB);
122 #ifdef HAVE_MUELU_EPETRA 123 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
125 RCP<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> >(vec);
126 if (tmpVec == Teuchos::null)
127 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
128 return tmpVec->getEpetra_MultiVector();
131 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
133 RCP<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> >(vec);
134 if (tmpVec == Teuchos::null)
135 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
136 return tmpVec->getEpetra_MultiVector();
139 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
141 const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> & tmpVec =
dynamic_cast<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> &
>(vec);
142 return *(tmpVec.getEpetra_MultiVector());
145 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
147 const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> & tmpVec =
dynamic_cast<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> &
>(vec);
148 return *(tmpVec.getEpetra_MultiVector());
151 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
153 RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<
const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
154 if (crsOp == Teuchos::null)
156 const RCP<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>>& tmp_ECrsMtx = rcp_dynamic_cast<
const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
157 if (tmp_ECrsMtx == Teuchos::null)
158 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
159 return tmp_ECrsMtx->getEpetra_CrsMatrix();
162 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
164 RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<
const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
165 if (crsOp == Teuchos::null)
167 const RCP<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<
const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
168 if (tmp_ECrsMtx == Teuchos::null)
169 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
170 return tmp_ECrsMtx->getEpetra_CrsMatrixNonConst();
173 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
176 const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>& crsOp =
dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(Op);
178 const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>& tmp_ECrsMtx =
dynamic_cast<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>&
>(*crsOp.getCrsMatrix());
179 return *tmp_ECrsMtx.getEpetra_CrsMatrix();
180 }
catch (std::bad_cast&) {
181 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
183 }
catch (std::bad_cast&) {
188 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
191 Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>& crsOp =
dynamic_cast<Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(Op);
193 Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>& tmp_ECrsMtx =
dynamic_cast<Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>&
>(*crsOp.getCrsMatrix());
194 return *tmp_ECrsMtx.getEpetra_CrsMatrixNonConst();
195 }
catch (std::bad_cast&) {
196 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
198 }
catch (std::bad_cast&) {
203 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
205 RCP<const Xpetra::EpetraMapT<GlobalOrdinal,Node> > xeMap = rcp_dynamic_cast<
const Xpetra::EpetraMapT<GlobalOrdinal,Node> >(rcpFromRef(map));
206 if (xeMap == Teuchos::null)
207 throw Exceptions::BadCast(
"Utilities::Map2EpetraMap : Cast from Xpetra::Map to Xpetra::EpetraMap failed");
208 return xeMap->getEpetra_Map();
212 #ifdef HAVE_MUELU_TPETRA 213 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
214 RCP<const Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
216 RCP<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(vec);
217 if (tmpVec == Teuchos::null)
218 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
219 return tmpVec->getTpetra_MultiVector();
222 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
224 RCP<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(vec);
225 if (tmpVec == Teuchos::null)
226 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
227 return tmpVec->getTpetra_MultiVector();
230 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
232 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec =
dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(vec);
233 return *(tmpVec.getTpetra_MultiVector());
236 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
238 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec =
dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(vec);
239 return tmpVec.getTpetra_MultiVector();
242 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
243 const Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>&
245 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec =
dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(vec);
246 return *(tmpVec.getTpetra_MultiVector());
249 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
252 RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<
const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
253 if (crsOp == Teuchos::null)
255 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<
const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
256 if (tmp_ECrsMtx == Teuchos::null)
257 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
258 return tmp_ECrsMtx->getTpetra_CrsMatrix();
261 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
263 RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<
const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
264 if (crsOp == Teuchos::null)
266 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<
const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
267 if (tmp_ECrsMtx == Teuchos::null)
268 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
269 return tmp_ECrsMtx->getTpetra_CrsMatrixNonConst();
272 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
275 const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>& crsOp =
dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(Op);
277 const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx =
dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(*crsOp.getCrsMatrix());
278 return *tmp_ECrsMtx.getTpetra_CrsMatrix();
279 }
catch (std::bad_cast&) {
280 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
282 }
catch (std::bad_cast&) {
287 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
290 Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>& crsOp =
dynamic_cast<Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(Op);
292 Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx =
dynamic_cast<Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(*crsOp.getCrsMatrix());
293 return *tmp_ECrsMtx.getTpetra_CrsMatrixNonConst();
294 }
catch (std::bad_cast&) {
295 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
297 }
catch (std::bad_cast&) {
302 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
304 RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<
const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
305 if (crsOp == Teuchos::null)
308 RCP<const Xpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsMat = crsOp->getCrsMatrix();
309 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_Crs = rcp_dynamic_cast<
const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
310 RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_BlockCrs;
311 if(!tmp_Crs.is_null()) {
312 return tmp_Crs->getTpetra_CrsMatrixNonConst();
315 tmp_BlockCrs= rcp_dynamic_cast<
const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
316 if (tmp_BlockCrs.is_null())
317 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
318 return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
322 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
324 RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<
const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
325 if (crsOp == Teuchos::null)
328 RCP<const Xpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsMat = crsOp->getCrsMatrix();
329 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_Crs = rcp_dynamic_cast<
const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
330 RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_BlockCrs;
331 if(!tmp_Crs.is_null()) {
332 return tmp_Crs->getTpetra_CrsMatrixNonConst();
335 tmp_BlockCrs= rcp_dynamic_cast<
const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
336 if (tmp_BlockCrs.is_null())
337 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
338 return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
343 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
345 const RCP<const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node>>& tmp_TMap = rcp_dynamic_cast<
const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node> >(rcpFromRef(map));
346 if (tmp_TMap == Teuchos::null)
347 throw Exceptions::BadCast(
"Utilities::Map2TpetraMap : Cast from Xpetra::Map to Xpetra::TpetraMap failed");
348 return tmp_TMap->getTpetra_Map();
352 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
355 bool doOptimizeStorage)
357 Scalar one = Teuchos::ScalarTraits<Scalar>::one();
358 Teuchos::ArrayRCP<Scalar> sv(scalingVector.size());
360 for (
int i = 0; i < scalingVector.size(); ++i)
361 sv[i] = one / scalingVector[i];
363 for (
int i = 0; i < scalingVector.size(); ++i)
364 sv[i] = scalingVector[i];
367 switch (Op.getRowMap()->lib()) {
368 case Xpetra::UseTpetra:
369 MyOldScaleMatrix_Tpetra(Op, sv, doFillComplete, doOptimizeStorage);
372 case Xpetra::UseEpetra:
373 MyOldScaleMatrix_Epetra(Op, sv, doFillComplete, doOptimizeStorage);
381 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
386 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
389 bool doOptimizeStorage)
391 #ifdef HAVE_MUELU_TPETRA 393 Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& tpOp = Op2NonConstTpetraCrs(Op);
395 const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rowMap = tpOp.getRowMap();
396 const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > domainMap = tpOp.getDomainMap();
397 const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rangeMap = tpOp.getRangeMap();
399 size_t maxRowSize = tpOp.getNodeMaxNumRowEntries();
400 if (maxRowSize == Teuchos::as<size_t>(-1))
403 std::vector<Scalar> scaledVals(maxRowSize);
404 if (tpOp.isFillComplete())
407 if (Op.isLocallyIndexed() ==
true) {
408 Teuchos::ArrayView<const LocalOrdinal> cols;
409 Teuchos::ArrayView<const Scalar> vals;
411 for (
size_t i = 0; i < rowMap->getNodeNumElements(); ++i) {
412 tpOp.getLocalRowView(i, cols, vals);
413 size_t nnz = tpOp.getNumEntriesInLocalRow(i);
414 if (nnz > maxRowSize) {
416 scaledVals.resize(maxRowSize);
418 for (
size_t j = 0; j < nnz; ++j)
419 scaledVals[j] = vals[j]*scalingVector[i];
422 Teuchos::ArrayView<const Scalar> valview(&scaledVals[0], nnz);
423 tpOp.replaceLocalValues(i, cols, valview);
428 Teuchos::ArrayView<const GlobalOrdinal> cols;
429 Teuchos::ArrayView<const Scalar> vals;
431 for (
size_t i = 0; i < rowMap->getNodeNumElements(); ++i) {
433 tpOp.getGlobalRowView(gid, cols, vals);
434 size_t nnz = tpOp.getNumEntriesInGlobalRow(gid);
435 if (nnz > maxRowSize) {
437 scaledVals.resize(maxRowSize);
440 for (
size_t j = 0; j < nnz; ++j)
441 scaledVals[j] = vals[j]*scalingVector[i];
444 Teuchos::ArrayView<const Scalar> valview(&scaledVals[0], nnz);
445 tpOp.replaceGlobalValues(gid, cols, valview);
450 if (doFillComplete) {
451 if (domainMap == Teuchos::null || rangeMap == Teuchos::null)
452 throw Exceptions::RuntimeError(
"In Utilities::Scaling: cannot fillComplete because the domain and/or range map hasn't been defined");
454 RCP<Teuchos::ParameterList> params = rcp(
new Teuchos::ParameterList());
455 params->set(
"Optimize Storage", doOptimizeStorage);
456 params->set(
"No Nonlocal Changes",
true);
457 Op.fillComplete(Op.getDomainMap(), Op.getRangeMap(), params);
467 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
468 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
470 Transpose (Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op,
bool ,
const std::string & label,
const Teuchos::RCP<Teuchos::ParameterList> ¶ms) {
471 #if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT) 472 std::string TorE =
"epetra";
474 std::string TorE =
"tpetra";
477 #if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT) 486 #ifdef HAVE_MUELU_TPETRA 487 if (TorE ==
"tpetra") {
491 RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > A;
492 Tpetra::RowMatrixTransposer<Scalar, LocalOrdinal, GlobalOrdinal, Node> transposer(rcpFromRef(tpetraOp),label);
495 using Teuchos::ParameterList;
497 RCP<ParameterList> transposeParams = params.is_null () ?
498 rcp (
new ParameterList) :
499 rcp (
new ParameterList (*params));
500 transposeParams->set (
"sort",
false);
501 A = transposer.createTranspose (transposeParams);
504 RCP<Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AA = rcp(
new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(A) );
505 RCP<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AAA = rcp_implicit_cast<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >(AA);
506 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AAAA = rcp(
new Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(AAA) );
507 if (!AAAA->isFillComplete())
508 AAAA->fillComplete(Op.getRangeMap(), Op.getDomainMap());
510 if (Op.IsView(
"stridedMaps"))
511 AAAA->CreateView(
"stridedMaps", Teuchos::rcpFromRef(Op),
true);
515 }
catch (std::exception& e) {
516 std::cout <<
"threw exception '" << e.what() <<
"'" << std::endl;
523 std::cout <<
"Utilities::Transpose() not implemented for Epetra" << std::endl;
524 return Teuchos::null;
529 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
530 RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
533 RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Xscalar;
534 #if defined(HAVE_XPETRA_TPETRA) && (defined(HAVE_TPETRA_INST_COMPLEX_DOUBLE) || defined(HAVE_TPETRA_INST_COMPLEX_FLOAT)) 535 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType real_type;
537 if ((
typeid(
Scalar).name() ==
typeid(std::complex<double>).name()) ||
538 (
typeid(
Scalar).name() ==
typeid(std::complex<float>).name())) {
539 Xscalar = Xpetra::MultiVectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(X->getMap(),X->getNumVectors());
540 size_t numVecs = X->getNumVectors();
541 for (
size_t j=0;j<numVecs;j++) {
542 Teuchos::ArrayRCP<const real_type> XVec = X->getData(j);
543 Teuchos::ArrayRCP<Scalar> XVecScalar = Xscalar->getDataNonConst(j);
544 for(
size_t i = 0; i < static_cast<size_t>(XVec.size()); ++i)
545 XVecScalar[i]=XVec[i];
549 Xscalar = rcp_dynamic_cast<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(X);
554 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
559 RCP<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,
LocalOrdinal,
GlobalOrdinal,
Node> > coordinates = Teuchos::null;
562 if(paramList.isParameter (
"Coordinates") ==
false)
565 #if defined(HAVE_MUELU_TPETRA) 571 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT) 572 typedef Tpetra::MultiVector<float, LocalOrdinal, GlobalOrdinal, Node> tfMV;
573 RCP<tfMV> floatCoords = Teuchos::null;
579 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_DOUBLE) 581 RCP<tdMV> doubleCoords = Teuchos::null;
582 if (paramList.isType<RCP<tdMV> >(
"Coordinates")) {
584 doubleCoords = paramList.get<RCP<tdMV> >(
"Coordinates");
585 paramList.remove(
"Coordinates");
587 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT) 588 else if (paramList.isType<RCP<tfMV> >(
"Coordinates")) {
590 floatCoords = paramList.get<RCP<tfMV> >(
"Coordinates");
591 paramList.remove(
"Coordinates");
592 doubleCoords = rcp(
new tdMV(floatCoords->getMap(), floatCoords->getNumVectors()));
597 if(doubleCoords != Teuchos::null) {
599 coordinates = Teuchos::rcp_dynamic_cast<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,
LocalOrdinal,
GlobalOrdinal,
Node> >(MueLu::TpetraMultiVector_To_XpetraMultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,
LocalOrdinal,
GlobalOrdinal,
Node>(doubleCoords));
600 TEUCHOS_TEST_FOR_EXCEPT(Teuchos::is_null(coordinates));
601 TEUCHOS_TEST_FOR_EXCEPT(doubleCoords->getNumVectors() != coordinates->getNumVectors());
606 throw Exceptions::RuntimeError(
"ExtractCoordinatesFromParameterList: The coordinates vector in parameter list is expected to be a Tpetra multivector with SC=double or float.");
608 #endif // endif HAVE_TPETRA 611 if(paramList.isType<decltype(coordinates)>(
"Coordinates")) {
612 coordinates = paramList.get<decltype(coordinates)>(
"Coordinates");
620 #define MUELU_UTILITIES_SHORT 621 #endif // MUELU_UTILITIES_DEF_HPP
Exception indicating invalid cast attempted.
static const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > Map2TpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
MueLu::DefaultLocalOrdinal LocalOrdinal
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 RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraCrs(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static const Epetra_Map & Map2EpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
static RCP< Epetra_CrsMatrix > Op2NonConstEpetraCrs(RCP< 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.
Namespace for MueLu classes and methods.
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec)
MueLu::DefaultScalar Scalar
MueLu::DefaultGlobalOrdinal GlobalOrdinal
static RCP< Epetra_MultiVector > MV2NonConstEpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec)
static RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > RealValuedToScalarMultiVector(RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::magnitudeType, LocalOrdinal, GlobalOrdinal, Node > > X)
static RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraCrs(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV2(Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
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 Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraRow(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
RCP< Xpetra::CrsMatrixWrap< SC, LO, GO, NO > > Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap(RCP< Epetra_CrsMatrix > &epAB)
void deep_copy(const DynRankView< DT, DP... > &dst, typename ViewTraits< DT, DP... >::const_value_type &value, typename std::enable_if< std::is_same< typename ViewTraits< DT, DP... >::specialize, void >::value >::type *=nullptr)
static void MyOldScaleMatrix_Epetra(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
static RCP< const Epetra_CrsMatrix > Op2EpetraCrs(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
Exception throws to report errors in the internal logical of the program.
static RCP< Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraRow(RCP< Xpetra::Matrix< 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)