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) 61 #include <EpetraExt_MatrixMatrix.h> 62 #include <EpetraExt_RowMatrixOut.h> 63 #include <EpetraExt_MultiVectorOut.h> 64 #include <EpetraExt_CrsMatrixIn.h> 65 #include <EpetraExt_MultiVectorIn.h> 66 #include <EpetraExt_BlockMapIn.h> 67 #include <Xpetra_EpetraUtils.hpp> 68 #include <Xpetra_EpetraMultiVector.hpp> 69 #include <EpetraExt_BlockMapOut.h> 72 #ifdef HAVE_MUELU_TPETRA 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_Import.hpp> 88 #include <Xpetra_ImportFactory.hpp> 89 #include <Xpetra_Map.hpp> 90 #include <Xpetra_MapFactory.hpp> 91 #include <Xpetra_Matrix.hpp> 92 #include <Xpetra_MatrixFactory.hpp> 93 #include <Xpetra_MultiVector.hpp> 94 #include <Xpetra_MultiVectorFactory.hpp> 95 #include <Xpetra_Operator.hpp> 96 #include <Xpetra_Vector.hpp> 97 #include <Xpetra_VectorFactory.hpp> 99 #include <Xpetra_MatrixMatrix.hpp> 102 #if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_ML) 103 #include <ml_operator.h> 104 #include <ml_epetra_utils.h> 109 #ifdef HAVE_MUELU_EPETRA 114 #ifdef HAVE_MUELU_EPETRA 116 template<
typename SC,
typename LO,
typename GO,
typename NO>
120 #ifdef HAVE_MUELU_EPETRA 121 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
123 RCP<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> >(vec);
124 if (tmpVec == Teuchos::null)
125 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
126 return tmpVec->getEpetra_MultiVector();
129 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
131 RCP<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> >(vec);
132 if (tmpVec == Teuchos::null)
133 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
134 return tmpVec->getEpetra_MultiVector();
137 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
139 const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> & tmpVec =
dynamic_cast<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> &
>(vec);
140 return *(tmpVec.getEpetra_MultiVector());
143 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
145 const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> & tmpVec =
dynamic_cast<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> &
>(vec);
146 return *(tmpVec.getEpetra_MultiVector());
149 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
151 RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<
const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
152 if (crsOp == Teuchos::null)
154 const RCP<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>>& tmp_ECrsMtx = rcp_dynamic_cast<
const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
155 if (tmp_ECrsMtx == Teuchos::null)
156 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
157 return tmp_ECrsMtx->getEpetra_CrsMatrix();
160 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
162 RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<
const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
163 if (crsOp == Teuchos::null)
165 const RCP<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<
const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
166 if (tmp_ECrsMtx == Teuchos::null)
167 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
168 return tmp_ECrsMtx->getEpetra_CrsMatrixNonConst();
171 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
174 const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>& crsOp =
dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(Op);
176 const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>& tmp_ECrsMtx =
dynamic_cast<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>&
>(*crsOp.getCrsMatrix());
177 return *tmp_ECrsMtx.getEpetra_CrsMatrix();
178 }
catch (std::bad_cast) {
179 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
181 }
catch (std::bad_cast) {
186 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
189 Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>& crsOp =
dynamic_cast<Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(Op);
191 Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>& tmp_ECrsMtx =
dynamic_cast<Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>&
>(*crsOp.getCrsMatrix());
192 return *tmp_ECrsMtx.getEpetra_CrsMatrixNonConst();
193 }
catch (std::bad_cast) {
194 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
196 }
catch (std::bad_cast) {
201 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
203 RCP<const Xpetra::EpetraMapT<GlobalOrdinal,Node> > xeMap = rcp_dynamic_cast<
const Xpetra::EpetraMapT<GlobalOrdinal,Node> >(rcpFromRef(map));
204 if (xeMap == Teuchos::null)
205 throw Exceptions::BadCast(
"Utilities::Map2EpetraMap : Cast from Xpetra::Map to Xpetra::EpetraMap failed");
206 return xeMap->getEpetra_Map();
210 #ifdef HAVE_MUELU_TPETRA 211 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
212 RCP<const Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
214 RCP<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(vec);
215 if (tmpVec == Teuchos::null)
216 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
217 return tmpVec->getTpetra_MultiVector();
220 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
222 RCP<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(vec);
223 if (tmpVec == Teuchos::null)
224 throw Exceptions::BadCast(
"Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
225 return tmpVec->getTpetra_MultiVector();
228 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
230 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec =
dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(vec);
231 return *(tmpVec.getTpetra_MultiVector());
234 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
236 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec =
dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(vec);
237 return tmpVec.getTpetra_MultiVector();
240 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
243 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec =
dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(vec);
244 return *(tmpVec.getTpetra_MultiVector());
247 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
250 RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<
const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
251 if (crsOp == Teuchos::null)
253 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<
const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
254 if (tmp_ECrsMtx == Teuchos::null)
255 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
256 return tmp_ECrsMtx->getTpetra_CrsMatrix();
259 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
261 RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<
const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
262 if (crsOp == Teuchos::null)
264 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<
const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
265 if (tmp_ECrsMtx == Teuchos::null)
266 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
267 return tmp_ECrsMtx->getTpetra_CrsMatrixNonConst();
270 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
273 const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>& crsOp =
dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(Op);
275 const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx =
dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(*crsOp.getCrsMatrix());
276 return *tmp_ECrsMtx.getTpetra_CrsMatrix();
277 }
catch (std::bad_cast) {
278 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
280 }
catch (std::bad_cast) {
285 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
288 Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>& crsOp =
dynamic_cast<Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(Op);
290 Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx =
dynamic_cast<Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&
>(*crsOp.getCrsMatrix());
291 return *tmp_ECrsMtx.getTpetra_CrsMatrixNonConst();
292 }
catch (std::bad_cast) {
293 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
295 }
catch (std::bad_cast) {
300 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
302 RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<
const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
303 if (crsOp == Teuchos::null)
306 RCP<const Xpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsMat = crsOp->getCrsMatrix();
307 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_Crs = rcp_dynamic_cast<
const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
308 RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_BlockCrs;
309 if(!tmp_Crs.is_null()) {
310 return tmp_Crs->getTpetra_CrsMatrixNonConst();
313 tmp_BlockCrs= rcp_dynamic_cast<
const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
314 if (tmp_BlockCrs.is_null())
315 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
316 return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
320 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
322 RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<
const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
323 if (crsOp == Teuchos::null)
326 RCP<const Xpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsMat = crsOp->getCrsMatrix();
327 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_Crs = rcp_dynamic_cast<
const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
328 RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_BlockCrs;
329 if(!tmp_Crs.is_null()) {
330 return tmp_Crs->getTpetra_CrsMatrixNonConst();
333 tmp_BlockCrs= rcp_dynamic_cast<
const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
334 if (tmp_BlockCrs.is_null())
335 throw Exceptions::BadCast(
"Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
336 return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
341 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
343 const RCP<const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node>>& tmp_TMap = rcp_dynamic_cast<
const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node> >(rcpFromRef(map));
344 if (tmp_TMap == Teuchos::null)
345 throw Exceptions::BadCast(
"Utilities::Map2TpetraMap : Cast from Xpetra::Map to Xpetra::TpetraMap failed");
346 return tmp_TMap->getTpetra_Map();
350 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
353 bool doOptimizeStorage)
355 Scalar one = Teuchos::ScalarTraits<Scalar>::one();
356 Teuchos::ArrayRCP<Scalar> sv(scalingVector.size());
358 for (
int i = 0; i < scalingVector.size(); ++i)
359 sv[i] = one / scalingVector[i];
361 for (
int i = 0; i < scalingVector.size(); ++i)
362 sv[i] = scalingVector[i];
365 switch (Op.getRowMap()->lib()) {
366 case Xpetra::UseTpetra:
367 MyOldScaleMatrix_Tpetra(Op, sv, doFillComplete, doOptimizeStorage);
370 case Xpetra::UseEpetra:
371 MyOldScaleMatrix_Epetra(Op, sv, doFillComplete, doOptimizeStorage);
380 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
385 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
388 bool doOptimizeStorage)
390 #ifdef HAVE_MUELU_TPETRA 394 const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rowMap = tpOp.
getRowMap();
395 const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > domainMap = tpOp.
getDomainMap();
396 const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rangeMap = tpOp.
getRangeMap();
399 if (maxRowSize == Teuchos::as<size_t>(-1))
402 std::vector<Scalar> scaledVals(maxRowSize);
406 if (Op.isLocallyIndexed() ==
true) {
407 Teuchos::ArrayView<const LocalOrdinal> cols;
408 Teuchos::ArrayView<const Scalar> vals;
410 for (
size_t i = 0; i < rowMap->getNodeNumElements(); ++i) {
413 if (nnz > maxRowSize) {
415 scaledVals.resize(maxRowSize);
417 for (
size_t j = 0; j < nnz; ++j)
418 scaledVals[j] = vals[j]*scalingVector[i];
421 Teuchos::ArrayView<const Scalar> valview(&scaledVals[0], nnz);
427 Teuchos::ArrayView<const GlobalOrdinal> cols;
428 Teuchos::ArrayView<const Scalar> vals;
430 for (
size_t i = 0; i < rowMap->getNodeNumElements(); ++i) {
431 GlobalOrdinal gid = rowMap->getGlobalElement(i);
434 if (nnz > maxRowSize) {
436 scaledVals.resize(maxRowSize);
439 for (
size_t j = 0; j < nnz; ++j)
440 scaledVals[j] = vals[j]*scalingVector[i];
443 Teuchos::ArrayView<const Scalar> valview(&scaledVals[0], nnz);
449 if (doFillComplete) {
450 if (domainMap == Teuchos::null || rangeMap == Teuchos::null)
451 throw Exceptions::RuntimeError(
"In Utilities::Scaling: cannot fillComplete because the domain and/or range map hasn't been defined");
453 RCP<Teuchos::ParameterList> params = rcp(
new Teuchos::ParameterList());
454 params->set(
"Optimize Storage", doOptimizeStorage);
455 params->set(
"No Nonlocal Changes",
true);
456 Op.fillComplete(Op.getDomainMap(), Op.getRangeMap(), params);
466 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
467 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
469 Transpose (Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op,
bool optimizeTranspose,
const std::string & label) {
470 #if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT) 471 std::string TorE =
"epetra";
473 std::string TorE =
"tpetra";
476 #if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT) 485 #ifdef HAVE_MUELU_TPETRA 486 if (TorE ==
"tpetra") {
490 RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > A;
494 RCP<Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AA = rcp(
new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(A) );
495 RCP<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AAA = rcp_implicit_cast<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >(AA);
496 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AAAA = rcp(
new Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(AAA) );
497 if (!AAAA->isFillComplete())
498 AAAA->fillComplete(Op.getRangeMap(), Op.getDomainMap());
500 if (Op.IsView(
"stridedMaps"))
501 AAAA->CreateView(
"stridedMaps", Teuchos::rcpFromRef(Op),
true);
505 }
catch (std::exception& e) {
506 std::cout <<
"threw exception '" << e.what() <<
"'" << std::endl;
513 std::cout <<
"Utilities::Transpose() not implemented for Epetra" << std::endl;
514 return Teuchos::null;
519 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
520 RCP<Xpetra::MultiVector<double,LocalOrdinal,GlobalOrdinal,Node> >
524 RCP<Xpetra::MultiVector<double,LocalOrdinal,GlobalOrdinal,Node> > coordinates = Teuchos::null;
527 if(paramList.isParameter (
"Coordinates") ==
false)
530 #if defined(HAVE_MUELU_TPETRA) 536 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT) 538 RCP<tfMV> floatCoords = Teuchos::null;
544 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_DOUBLE) 546 RCP<tdMV> doubleCoords = Teuchos::null;
547 if (paramList.isType<RCP<tdMV> >(
"Coordinates")) {
549 doubleCoords = paramList.get<RCP<tdMV> >(
"Coordinates");
550 paramList.remove(
"Coordinates");
552 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT) 553 else if (paramList.isType<RCP<tfMV> >(
"Coordinates")) {
555 floatCoords = paramList.get<RCP<tfMV> >(
"Coordinates");
556 paramList.remove(
"Coordinates");
557 doubleCoords = rcp(
new tdMV(floatCoords->getMap(), floatCoords->getNumVectors()));
562 if(doubleCoords != Teuchos::null) {
564 coordinates = Teuchos::rcp_dynamic_cast<Xpetra::MultiVector<double,LocalOrdinal,GlobalOrdinal,Node> >(MueLu::TpetraMultiVector_To_XpetraMultiVector<double,LocalOrdinal,GlobalOrdinal,Node>(doubleCoords));
565 TEUCHOS_TEST_FOR_EXCEPT(Teuchos::is_null(coordinates));
566 TEUCHOS_TEST_FOR_EXCEPT(doubleCoords->getNumVectors() != coordinates->getNumVectors());
571 throw Exceptions::RuntimeError(
"ExtractCoordinatesFromParameterList: The coordinates vector in parameter list is expected to be a Tpetra multivector with SC=double or float.");
573 #endif // endif HAVE_TPETRA 580 #define MUELU_UTILITIES_SHORT 581 #endif // MUELU_UTILITIES_DEF_HPP LocalOrdinal replaceGlobalValues(const GlobalOrdinal globalRow, const typename UnmanagedView< GlobalIndicesViewType >::type &inputInds, const typename UnmanagedView< ImplScalarViewType >::type &inputVals) const
void getGlobalRowView(GlobalOrdinal GlobalRow, Teuchos::ArrayView< const GlobalOrdinal > &indices, Teuchos::ArrayView< const Scalar > &values) const
Exception indicating invalid cast attempted.
static const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > Map2TpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
size_t getNodeMaxNumRowEntries() const
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.
LocalOrdinal replaceLocalValues(const LocalOrdinal localRow, const typename UnmanagedView< LocalIndicesViewType >::type &inputInds, const typename UnmanagedView< ImplScalarViewType >::type &inputVals) const
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< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Transpose(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, bool optimizeTranspose=false, const std::string &label=std::string())
static RCP< Epetra_CrsMatrix > Op2NonConstEpetraCrs(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
void resumeFill(const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
Teuchos::RCP< const map_type > getRowMap() const
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.
size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
void deep_copy(MultiVector< DS, DL, DG, DN, dstClassic > &dst, const MultiVector< SS, SL, SG, SN, srcClassic > &src)
Namespace for MueLu classes and methods.
static RCP< Xpetra::MultiVector< double, LocalOrdinal, GlobalOrdinal, Node > > ExtractCoordinatesFromParameterList(ParameterList ¶mList)
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)
bool isFillComplete() const
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< const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraRow(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
Teuchos::RCP< crs_matrix_type > createTranspose()
void getLocalRowView(LocalOrdinal LocalRow, Teuchos::ArrayView< const LocalOrdinal > &indices, Teuchos::ArrayView< const Scalar > &values) const
RCP< Xpetra::CrsMatrixWrap< SC, LO, GO, NO > > Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap(RCP< Epetra_CrsMatrix > &epAB)
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
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)
Teuchos::RCP< const map_type > getDomainMap() const
Exception throws to report errors in the internal logical of the program.
Teuchos::RCP< const map_type > getRangeMap() const
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)