MueLu  Version of the Day
MueLu_Utilities_decl.hpp
Go to the documentation of this file.
1 // @HEADER
2 //
3 // ***********************************************************************
4 //
5 // MueLu: A package for multigrid based preconditioning
6 // Copyright 2012 Sandia Corporation
7 //
8 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9 // the U.S. Government retains certain rights in this software.
10 //
11 // Redistribution and use in source and binary forms, with or without
12 // modification, are permitted provided that the following conditions are
13 // met:
14 //
15 // 1. Redistributions of source code must retain the above copyright
16 // notice, this list of conditions and the following disclaimer.
17 //
18 // 2. Redistributions in binary form must reproduce the above copyright
19 // notice, this list of conditions and the following disclaimer in the
20 // documentation and/or other materials provided with the distribution.
21 //
22 // 3. Neither the name of the Corporation nor the names of the
23 // contributors may be used to endorse or promote products derived from
24 // this software without specific prior written permission.
25 //
26 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37 //
38 // Questions? Contact
39 // Jonathan Hu (jhu@sandia.gov)
40 // Andrey Prokopenko (aprokop@sandia.gov)
41 // Ray Tuminaro (rstumin@sandia.gov)
42 //
43 // ***********************************************************************
44 //
45 // @HEADER
46 #ifndef MUELU_UTILITIES_DECL_HPP
47 #define MUELU_UTILITIES_DECL_HPP
48 
49 #include <unistd.h> //necessary for "sleep" function in debugging methods
50 #include <string>
51 
52 #include "MueLu_ConfigDefs.hpp"
53 
54 #include <Teuchos_DefaultComm.hpp>
55 #include <Teuchos_ScalarTraits.hpp>
56 #include <Teuchos_ParameterList.hpp>
57 
58 #ifdef HAVE_MUELU_TPETRA
60 #endif
62 #include <Xpetra_CrsMatrix_fwd.hpp>
64 #include <Xpetra_Map_fwd.hpp>
66 #include <Xpetra_Matrix_fwd.hpp>
70 #include <Xpetra_Operator_fwd.hpp>
71 #include <Xpetra_Vector_fwd.hpp>
73 #include <Xpetra_ExportFactory.hpp>
74 
75 #include <Xpetra_Import.hpp>
76 #include <Xpetra_ImportFactory.hpp>
77 #include <Xpetra_MatrixMatrix.hpp>
78 
79 #ifdef HAVE_MUELU_EPETRA
81 
82 // needed because of inlined function
83 //TODO: remove inline function?
85 #include <Xpetra_CrsMatrixWrap.hpp>
86 
87 #endif
88 
89 #include "MueLu_Exceptions.hpp"
90 
91 #ifdef HAVE_MUELU_EPETRAEXT
92 class Epetra_CrsMatrix;
93 class Epetra_MultiVector;
94 class Epetra_Vector;
95 #include "EpetraExt_Transpose_RowMatrix.h"
96 #endif
97 
98 #ifdef HAVE_MUELU_TPETRA
99 #include <Tpetra_CrsMatrix.hpp>
100 #include <Tpetra_RowMatrixTransposer.hpp>
101 #include <Tpetra_Map.hpp>
102 #include <Tpetra_MultiVector.hpp>
105 #endif
106 
107 #include <MueLu_UtilitiesBase.hpp>
108 
109 
110 namespace MueLu {
111 
112 #ifdef HAVE_MUELU_EPETRA
113  //defined after Utilities class
114  template<typename SC,typename LO,typename GO,typename NO>
115  RCP<Xpetra::CrsMatrixWrap<SC,LO,GO,NO> >
116  Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap(RCP<Epetra_CrsMatrix> &epAB);
117 
118  template<typename SC,typename LO,typename GO,typename NO>
119  RCP<Xpetra::Matrix<SC, LO, GO, NO> >
120  EpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Epetra_CrsMatrix>& A);
121 
122  template<typename SC,typename LO,typename GO,typename NO>
123  RCP<Xpetra::MultiVector<SC, LO, GO, NO> >
124  EpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Epetra_MultiVector>& V);
125 #endif
126 
127 #ifdef HAVE_MUELU_TPETRA
128  template<typename SC,typename LO,typename GO,typename NO>
129  RCP<Xpetra::Matrix<SC, LO, GO, NO> >
130  TpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Tpetra::CrsMatrix<SC, LO, GO, NO> >& Atpetra);
131 
132 
133  template<typename SC,typename LO,typename GO,typename NO>
134  RCP<Xpetra::MultiVector<SC, LO, GO, NO> >
136 #endif
137 
145  template <class Scalar,
146  class LocalOrdinal = int,
147  class GlobalOrdinal = LocalOrdinal,
149  class Utilities : public UtilitiesBase<Scalar, LocalOrdinal, GlobalOrdinal, Node> {
150 #undef MUELU_UTILITIES_SHORT
151  //#include "MueLu_UseShortNames.hpp"
152 
153  public:
154  typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType Magnitude;
155 
156 #ifdef HAVE_MUELU_EPETRA
157  // @{
159  static RCP<const Epetra_MultiVector> MV2EpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > const vec);
160  static RCP< Epetra_MultiVector> MV2NonConstEpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > vec);
161 
162  static const Epetra_MultiVector& MV2EpetraMV(const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
164 
165  static RCP<const Epetra_CrsMatrix> Op2EpetraCrs(RCP<const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
167 
168  static const Epetra_CrsMatrix& Op2EpetraCrs(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
170 
171  static const Epetra_Map& Map2EpetraMap(const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node>& map);
172  // @}
173 #endif
174 
175 #ifdef HAVE_MUELU_TPETRA
176  static RCP<const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2TpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > const vec);
178  static RCP< Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2NonConstTpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > vec);
179  static RCP< Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2NonConstTpetraMV2(Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
180 
183 
184  static RCP<const Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2TpetraCrs(RCP<const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
185  static RCP< Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2NonConstTpetraCrs(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
186 
189 
190  static RCP<const Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2TpetraRow(RCP<const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
191  static RCP< Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2NonConstTpetraRow(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
192 
193 
194  static const RCP<const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > Map2TpetraMap(const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node>& map);
195 #endif
196 
197  static RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Crs2Op(RCP<Xpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op) { return UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Crs2Op(Op); }
199  static RCP<Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > GetMatrixDiagonalInverse(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, Magnitude tol = Teuchos::ScalarTraits<Scalar>::eps()*100) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetMatrixDiagonalInverse(A,tol); }
201  static Teuchos::RCP<Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > GetLumpedMatrixDiagonal(Teuchos::RCP<const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > A) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetLumpedMatrixDiagonal(A); }
203  static Teuchos::RCP<Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > GetInverse(Teuchos::RCP<const Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > v, Magnitude tol = Teuchos::ScalarTraits<Scalar>::eps()*100, Scalar tolReplacement = Teuchos::ScalarTraits<Scalar>::zero()) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetInverse(v,tol,tolReplacement); }
207  static RCP<Teuchos::FancyOStream> MakeFancy(std::ostream& os) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::MakeFancy(os); }
208  static typename Teuchos::ScalarTraits<Scalar>::magnitudeType Distance2(const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& v, LocalOrdinal i0, LocalOrdinal i1) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Distance2(v,i0,i1); }
209  static Teuchos::ArrayRCP<const bool> DetectDirichletRows(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, const Magnitude& tol = Teuchos::ScalarTraits<Scalar>::magnitude(0.)) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::DetectDirichletRows(A,tol); }
211 
212  static Scalar PowerMethod(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, bool scaleByDiag = true,
213  LocalOrdinal niters = 10, Magnitude tolerance = 1e-2, bool verbose = false, unsigned int seed = 123) {
214  return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::PowerMethod(A,scaleByDiag,niters,tolerance,verbose,seed);
215  }
216 
219  }
220 
221  static void MyOldScaleMatrix(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Teuchos::ArrayRCP<const Scalar>& scalingVector, bool doInverse = true,
222  bool doFillComplete = true, bool doOptimizeStorage = true);
223 
224  static void MyOldScaleMatrix_Epetra(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Teuchos::ArrayRCP<Scalar>& scalingVector,
225  bool doFillComplete, bool doOptimizeStorage);
226  static void MyOldScaleMatrix_Tpetra(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Teuchos::ArrayRCP<Scalar>& scalingVector,
227  bool doFillComplete, bool doOptimizeStorage);
228 
229  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());
230 
231  static RCP<Xpetra::MultiVector<double,LocalOrdinal,GlobalOrdinal,Node> > ExtractCoordinatesFromParameterList(ParameterList& paramList);
232  }; // class Utilities
233 
235 
236 #ifdef HAVE_MUELU_EPETRA
237 
246  template <>
247  class Utilities<double,int,int,Xpetra::EpetraNode> : public UtilitiesBase<double,int,int,Xpetra::EpetraNode> {
248  public:
249  typedef double Scalar;
250  typedef int LocalOrdinal;
251  typedef int GlobalOrdinal;
253  typedef Teuchos::ScalarTraits<Scalar>::magnitudeType Magnitude;
254 
255  private:
262 #ifdef HAVE_MUELU_EPETRA
266 #endif
267  public:
268 
269 #ifdef HAVE_MUELU_EPETRA
270  // @{
272  static RCP<const Epetra_MultiVector> MV2EpetraMV(RCP<MultiVector> const vec) {
273  RCP<const EpetraMultiVector > tmpVec = rcp_dynamic_cast<EpetraMultiVector>(vec);
274  if (tmpVec == Teuchos::null)
275  throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
276  return tmpVec->getEpetra_MultiVector();
277  }
278  static RCP< Epetra_MultiVector> MV2NonConstEpetraMV(RCP<MultiVector> vec) {
279  RCP<const EpetraMultiVector> tmpVec = rcp_dynamic_cast<EpetraMultiVector>(vec);
280  if (tmpVec == Teuchos::null)
281  throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
282  return tmpVec->getEpetra_MultiVector();
283  }
284 
285  static const Epetra_MultiVector& MV2EpetraMV(const MultiVector& vec) {
286  const EpetraMultiVector& tmpVec = dynamic_cast<const EpetraMultiVector&>(vec);
287  return *(tmpVec.getEpetra_MultiVector());
288  }
289  static Epetra_MultiVector& MV2NonConstEpetraMV(MultiVector& vec) {
290  const EpetraMultiVector& tmpVec = dynamic_cast<const EpetraMultiVector&>(vec);
291  return *(tmpVec.getEpetra_MultiVector());
292  }
293 
294  static RCP<const Epetra_CrsMatrix> Op2EpetraCrs(RCP<const Matrix> Op) {
295  RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
296  if (crsOp == Teuchos::null)
297  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
298  const RCP<const EpetraCrsMatrix>& tmp_ECrsMtx = rcp_dynamic_cast<const EpetraCrsMatrix>(crsOp->getCrsMatrix());
299  if (tmp_ECrsMtx == Teuchos::null)
300  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
301  return tmp_ECrsMtx->getEpetra_CrsMatrix();
302  }
303  static RCP< Epetra_CrsMatrix> Op2NonConstEpetraCrs(RCP<Matrix> Op) {
304  RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
305  if (crsOp == Teuchos::null)
306  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
307  const RCP<const EpetraCrsMatrix> &tmp_ECrsMtx = rcp_dynamic_cast<const EpetraCrsMatrix>(crsOp->getCrsMatrix());
308  if (tmp_ECrsMtx == Teuchos::null)
309  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
310  return tmp_ECrsMtx->getEpetra_CrsMatrixNonConst();
311  }
312 
313  static const Epetra_CrsMatrix& Op2EpetraCrs(const Matrix& Op) {
314  try {
315  const CrsMatrixWrap& crsOp = dynamic_cast<const CrsMatrixWrap&>(Op);
316  try {
317  const EpetraCrsMatrix& tmp_ECrsMtx = dynamic_cast<const EpetraCrsMatrix&>(*crsOp.getCrsMatrix());
318  return *tmp_ECrsMtx.getEpetra_CrsMatrix();
319  } catch (std::bad_cast) {
320  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
321  }
322  } catch (std::bad_cast) {
323  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
324  }
325  }
326  static Epetra_CrsMatrix& Op2NonConstEpetraCrs(Matrix& Op) {
327  try {
328  CrsMatrixWrap& crsOp = dynamic_cast<CrsMatrixWrap&>(Op);
329  try {
330  EpetraCrsMatrix& tmp_ECrsMtx = dynamic_cast<EpetraCrsMatrix&>(*crsOp.getCrsMatrix());
331  return *tmp_ECrsMtx.getEpetra_CrsMatrixNonConst();
332  } catch (std::bad_cast) {
333  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
334  }
335  } catch (std::bad_cast) {
336  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
337  }
338  }
339 
340  static const Epetra_Map& Map2EpetraMap(const Map& map) {
341  RCP<const EpetraMap> xeMap = rcp_dynamic_cast<const EpetraMap>(rcpFromRef(map));
342  if (xeMap == Teuchos::null)
343  throw Exceptions::BadCast("Utilities::Map2EpetraMap : Cast from Xpetra::Map to Xpetra::EpetraMap failed");
344  return xeMap->getEpetra_Map();
345  }
346  // @}
347 #endif
348 
349 #ifdef HAVE_MUELU_TPETRA
350  static RCP<const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2TpetraMV(RCP<MultiVector> const vec) {
352 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
353  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
354  throw Exceptions::RuntimeError("MV2TpetraMV: Tpetra has not been compiled with support for LO=GO=int.");
355 #else
356  RCP<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(vec);
357  if (tmpVec == Teuchos::null)
358  throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
359  return tmpVec->getTpetra_MultiVector();
360 #endif
361  }
362  static RCP< Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2NonConstTpetraMV(RCP<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))))
365  throw Exceptions::RuntimeError("MV2NonConstTpetraMV: Tpetra has not been compiled with support for LO=GO=int.");
366 #else
367  RCP<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(vec);
368  if (tmpVec == Teuchos::null)
369  throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
370  return tmpVec->getTpetra_MultiVector();
371 #endif
372 
373  }
374  static RCP< Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2NonConstTpetraMV2(MultiVector& vec) {
375 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
376  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
377  throw Exceptions::RuntimeError("MV2NonConstTpetraMV2: Tpetra has not been compiled with support for LO=GO=int.");
378 #else
380  return tmpVec.getTpetra_MultiVector();
381 #endif
382  }
383 
385 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
386  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
387  throw Exceptions::RuntimeError("MV2TpetraMV: Tpetra has not been compiled with support for LO=GO=int.");
388 #else
390  return *(tmpVec.getTpetra_MultiVector());
391 #endif
392  }
394 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
395  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
396  throw Exceptions::RuntimeError("MV2NonConstTpetraMV: Tpetra has not been compiled with support for LO=GO=int.");
397 #else
399  return *(tmpVec.getTpetra_MultiVector());
400 #endif
401  }
402 
403  static RCP<const Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2TpetraCrs(RCP<const Matrix> Op) {
404 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
405  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
406  throw Exceptions::RuntimeError("Op2TpetraCrs: Tpetra has not been compiled with support for LO=GO=int.");
407 #else
408  // Get the underlying Tpetra Mtx
409  RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
410  if (crsOp == Teuchos::null)
411  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
412  const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
413  if (tmp_ECrsMtx == Teuchos::null)
414  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
415  return tmp_ECrsMtx->getTpetra_CrsMatrix();
416 #endif
417  }
418  static RCP< Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2NonConstTpetraCrs(RCP<Matrix> Op){
419 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
420  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
421  throw Exceptions::RuntimeError("Op2NonConstTpetraCrs: Tpetra has not been compiled with support for LO=GO=int.");
422 #else
423  RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
424  if (crsOp == Teuchos::null)
425  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
426  const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
427  if (tmp_ECrsMtx == Teuchos::null)
428  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
429  return tmp_ECrsMtx->getTpetra_CrsMatrixNonConst();
430 #endif
431  };
432 
434 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
435  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
436  throw Exceptions::RuntimeError("Op2TpetraCrs: Tpetra has not been compiled with support for LO=GO=int.");
437 #else
438  try {
439  const CrsMatrixWrap& crsOp = dynamic_cast<const CrsMatrixWrap&>(Op);
440  try {
442  return *tmp_ECrsMtx.getTpetra_CrsMatrix();
443  } catch (std::bad_cast) {
444  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
445  }
446  } catch (std::bad_cast) {
447  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
448  }
449 #endif
450  }
452 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
453  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
454  throw Exceptions::RuntimeError("Op2NonConstTpetraCrs: Tpetra has not been compiled with support for LO=GO=int.");
455 #else
456  try {
457  CrsMatrixWrap& crsOp = dynamic_cast<CrsMatrixWrap&>(Op);
458  try {
460  return *tmp_ECrsMtx.getTpetra_CrsMatrixNonConst();
461  } catch (std::bad_cast) {
462  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
463  }
464  } catch (std::bad_cast) {
465  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
466  }
467 #endif
468  }
469 
470  static RCP<const Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2TpetraRow(RCP<const Matrix> Op) {
471 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
472  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
473  throw Exceptions::RuntimeError("Op2TpetraRow: Tpetra has not been compiled with support for LO=GO=int.");
474 #else
475  RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
476  if (crsOp == Teuchos::null)
477  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
478 
479  RCP<const CrsMatrix> crsMat = crsOp->getCrsMatrix();
480  const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_Crs = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
481  RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_BlockCrs;
482  if(!tmp_Crs.is_null()) {
483  return tmp_Crs->getTpetra_CrsMatrixNonConst();
484  }
485  else {
486  tmp_BlockCrs= rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
487  if (tmp_BlockCrs.is_null())
488  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
489  return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
490  }
491 #endif
492  }
493  static RCP< Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2NonConstTpetraRow(RCP<Matrix> Op) {
494 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
495  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
496  throw Exceptions::RuntimeError("Op2NonConstTpetraRow: Tpetra has not been compiled with support for LO=GO=int.");
497 #else
498  RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
499  if (crsOp == Teuchos::null)
500  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
501 
502  RCP<const CrsMatrix> crsMat = crsOp->getCrsMatrix();
503  const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_Crs = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
504  RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_BlockCrs;
505  if(!tmp_Crs.is_null()) {
506  return tmp_Crs->getTpetra_CrsMatrixNonConst();
507  }
508  else {
509  tmp_BlockCrs= rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
510  if (tmp_BlockCrs.is_null())
511  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
512  return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
513  }
514 #endif
515  };
516 
517 
518  static const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > Map2TpetraMap(const Map& map) {
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))))
521  throw Exceptions::RuntimeError("Map2TpetraMap: Tpetra has not been compiled with support for LO=GO=int.");
522 #else
523  const RCP<const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node>>& tmp_TMap = rcp_dynamic_cast<const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node> >(rcpFromRef(map));
524  if (tmp_TMap == Teuchos::null)
525  throw Exceptions::BadCast("Utilities::Map2TpetraMap : Cast from Xpetra::Map to Xpetra::TpetraMap failed");
526  return tmp_TMap->getTpetra_Map();
527 #endif
528  };
529 #endif
530 
531  static RCP<Matrix> Crs2Op(RCP<CrsMatrix> Op) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Crs2Op(Op); }
533  static RCP<Vector> GetMatrixDiagonalInverse(const Matrix& A, Magnitude tol = Teuchos::ScalarTraits<Scalar>::eps()*100) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetMatrixDiagonalInverse(A,tol); }
535  static Teuchos::RCP<Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > GetLumpedMatrixDiagonal(Teuchos::RCP<const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > A) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetLumpedMatrixDiagonal(A); }
537  static RCP<Vector> GetInverse(Teuchos::RCP<const Vector> v, Magnitude tol = Teuchos::ScalarTraits<Scalar>::eps()*100, Scalar tolReplacement = Teuchos::ScalarTraits<Scalar>::zero()) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetInverse(v,tol,tolReplacement); }
541  static RCP<Teuchos::FancyOStream> MakeFancy(std::ostream& os) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::MakeFancy(os); }
542  static Teuchos::ScalarTraits<Scalar>::magnitudeType Distance2(const MultiVector& v, LocalOrdinal i0, LocalOrdinal i1) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Distance2(v,i0,i1); }
543  static Teuchos::ArrayRCP<const bool> DetectDirichletRows(const Matrix& A, const Magnitude& tol = Teuchos::ScalarTraits<Scalar>::zero()) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::DetectDirichletRows(A,tol); }
545 
546  static Scalar PowerMethod(const Matrix& A, bool scaleByDiag = true,
547  LocalOrdinal niters = 10, Magnitude tolerance = 1e-2, bool verbose = false, unsigned int seed = 123) {
548  return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::PowerMethod(A,scaleByDiag,niters,tolerance,verbose,seed);
549  }
550 
553  }
554 
555  static void MyOldScaleMatrix(Matrix& Op, const Teuchos::ArrayRCP<const Scalar>& scalingVector, bool doInverse = true,
556  bool doFillComplete = true, bool doOptimizeStorage = true) {
557  Scalar one = Teuchos::ScalarTraits<Scalar>::one();
558  Teuchos::ArrayRCP<Scalar> sv(scalingVector.size());
559  if (doInverse) {
560  for (int i = 0; i < scalingVector.size(); ++i)
561  sv[i] = one / scalingVector[i];
562  } else {
563  for (int i = 0; i < scalingVector.size(); ++i)
564  sv[i] = scalingVector[i];
565  }
566 
567  switch (Op.getRowMap()->lib()) {
568  case Xpetra::UseTpetra:
569  MyOldScaleMatrix_Tpetra(Op, sv, doFillComplete, doOptimizeStorage);
570  break;
571 
572  case Xpetra::UseEpetra:
573  MyOldScaleMatrix_Epetra(Op, sv, doFillComplete, doOptimizeStorage);
574  break;
575 
576  default:
577  throw Exceptions::RuntimeError("Only Epetra and Tpetra matrices can be scaled.");
578  break;
579  }
580  }
581 
582  // TODO This is the <double,int,int> specialization
583  static void MyOldScaleMatrix_Tpetra(Matrix& Op, const Teuchos::ArrayRCP<Scalar>& scalingVector,
584  bool doFillComplete, bool doOptimizeStorage) {
585 #ifdef HAVE_MUELU_TPETRA
586 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
587  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
588  throw Exceptions::RuntimeError("Matrix scaling is not possible because Tpetra has not been compiled with support for LO=GO=int.");
589 #else
590  try {
592 
593  const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rowMap = tpOp.getRowMap();
594  const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > domainMap = tpOp.getDomainMap();
595  const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rangeMap = tpOp.getRangeMap();
596 
597  size_t maxRowSize = tpOp.getNodeMaxNumRowEntries();
598  if (maxRowSize == Teuchos::as<size_t>(-1)) // hasn't been determined yet
599  maxRowSize = 20;
600 
601  std::vector<Scalar> scaledVals(maxRowSize);
602  if (tpOp.isFillComplete())
603  tpOp.resumeFill();
604 
605  if (Op.isLocallyIndexed() == true) {
606  Teuchos::ArrayView<const LocalOrdinal> cols;
607  Teuchos::ArrayView<const Scalar> vals;
608 
609  for (size_t i = 0; i < rowMap->getNodeNumElements(); ++i) {
610  tpOp.getLocalRowView(i, cols, vals);
611  size_t nnz = tpOp.getNumEntriesInLocalRow(i);
612  if (nnz > maxRowSize) {
613  maxRowSize = nnz;
614  scaledVals.resize(maxRowSize);
615  }
616  for (size_t j = 0; j < nnz; ++j)
617  scaledVals[j] = vals[j]*scalingVector[i];
618 
619  if (nnz > 0) {
620  Teuchos::ArrayView<const Scalar> valview(&scaledVals[0], nnz);
621  tpOp.replaceLocalValues(i, cols, valview);
622  }
623  } //for (size_t i=0; ...
624 
625  } else {
626  Teuchos::ArrayView<const GlobalOrdinal> cols;
627  Teuchos::ArrayView<const Scalar> vals;
628 
629  for (size_t i = 0; i < rowMap->getNodeNumElements(); ++i) {
630  GlobalOrdinal gid = rowMap->getGlobalElement(i);
631  tpOp.getGlobalRowView(gid, cols, vals);
632  size_t nnz = tpOp.getNumEntriesInGlobalRow(gid);
633  if (nnz > maxRowSize) {
634  maxRowSize = nnz;
635  scaledVals.resize(maxRowSize);
636  }
637  // FIXME FIXME FIXME FIXME FIXME FIXME
638  for (size_t j = 0; j < nnz; ++j)
639  scaledVals[j] = vals[j]*scalingVector[i]; //FIXME i or gid?
640 
641  if (nnz > 0) {
642  Teuchos::ArrayView<const Scalar> valview(&scaledVals[0], nnz);
643  tpOp.replaceGlobalValues(gid, cols, valview);
644  }
645  } //for (size_t i=0; ...
646  }
647 
648  if (doFillComplete) {
649  if (domainMap == Teuchos::null || rangeMap == Teuchos::null)
650  throw Exceptions::RuntimeError("In Utilities::Scaling: cannot fillComplete because the domain and/or range map hasn't been defined");
651 
652  RCP<Teuchos::ParameterList> params = rcp(new Teuchos::ParameterList());
653  params->set("Optimize Storage", doOptimizeStorage);
654  params->set("No Nonlocal Changes", true);
655  Op.fillComplete(Op.getDomainMap(), Op.getRangeMap(), params);
656  }
657  } catch(...) {
658  throw Exceptions::RuntimeError("Only Tpetra::CrsMatrix types can be scaled (Err.1)");
659  }
660 #endif
661 #else
662  throw Exceptions::RuntimeError("Matrix scaling is not possible because Tpetra has not been enabled.");
663 #endif
664  }
665 
666  static void MyOldScaleMatrix_Epetra (Matrix& Op, const Teuchos::ArrayRCP<Scalar>& scalingVector, bool doFillComplete, bool doOptimizeStorage) {
667 #ifdef HAVE_MUELU_EPETRA
668  try {
669  //const Epetra_CrsMatrix& epOp = Utilities<double,int,int>::Op2NonConstEpetraCrs(Op);
670  const Epetra_CrsMatrix& epOp = Op2NonConstEpetraCrs(Op);
671 
672  Epetra_Map const &rowMap = epOp.RowMap();
673  int nnz;
674  double *vals;
675  int *cols;
676 
677  for (int i = 0; i < rowMap.NumMyElements(); ++i) {
678  epOp.ExtractMyRowView(i, nnz, vals, cols);
679  for (int j = 0; j < nnz; ++j)
680  vals[j] *= scalingVector[i];
681  }
682 
683  } catch (...){
684  throw Exceptions::RuntimeError("Only Epetra_CrsMatrix types can be scaled");
685  }
686 #else
687  throw Exceptions::RuntimeError("Matrix scaling is not possible because Epetra has not been enabled.");
688 #endif // HAVE_MUELU_EPETRA
689  }
690 
696  static RCP<Matrix> Transpose(Matrix& Op, bool optimizeTranspose = false,const std::string & label = std::string()) {
697  switch (Op.getRowMap()->lib()) {
698  case Xpetra::UseTpetra: {
699 #ifdef HAVE_MUELU_TPETRA
700 #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
701  (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
702  throw Exceptions::RuntimeError("Utilities::Transpose: Tpetra is not compiled with LO=GO=int. Add TPETRA_INST_INT_INT:BOOL=ON to your configuration!");
703 #else
704  try {
706 
707  // Compute the transpose A of the Tpetra matrix tpetraOp.
708  RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > A;
709  Tpetra::RowMatrixTransposer<Scalar, LocalOrdinal, GlobalOrdinal, Node> transposer(rcpFromRef(tpetraOp),label);
710  A = transposer.createTranspose();
711  RCP<Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AA = rcp(new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(A));
712  RCP<CrsMatrix> AAA = rcp_implicit_cast<CrsMatrix>(AA);
713  RCP<Matrix> AAAA = rcp( new CrsMatrixWrap(AAA));
714 
715  if (Op.IsView("stridedMaps"))
716  AAAA->CreateView("stridedMaps", Teuchos::rcpFromRef(Op), true/*doTranspose*/);
717 
718  return AAAA;
719  }
720  catch (std::exception& e) {
721  std::cout << "threw exception '" << e.what() << "'" << std::endl;
722  throw Exceptions::RuntimeError("Utilities::Transpose failed, perhaps because matrix is not a Crs matrix");
723  }
724 #endif
725 #else
726  throw Exceptions::RuntimeError("Utilities::Transpose: Tpetra is not compiled!");
727 #endif
728  break;
729  }
730  case Xpetra::UseEpetra:
731  {
732 #if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT)
733  Teuchos::TimeMonitor tm(*Teuchos::TimeMonitor::getNewTimer("ZZ Entire Transpose"));
734  // Epetra case
736  EpetraExt::RowMatrix_Transpose transposer;
737  Epetra_CrsMatrix * A = dynamic_cast<Epetra_CrsMatrix*>(&transposer(epetraOp));
738  transposer.ReleaseTranspose(); // So we can keep A in Muelu...
739 
740  RCP<Epetra_CrsMatrix> rcpA(A);
741  RCP<EpetraCrsMatrix> AA = rcp(new EpetraCrsMatrix(rcpA));
742  RCP<CrsMatrix> AAA = rcp_implicit_cast<CrsMatrix>(AA);
743  RCP<Matrix> AAAA = rcp( new CrsMatrixWrap(AAA));
744  AAAA->fillComplete(Op.getRangeMap(), Op.getDomainMap());
745 
746  if (Op.IsView("stridedMaps"))
747  AAAA->CreateView("stridedMaps", Teuchos::rcpFromRef(Op), true/*doTranspose*/);
748 
749  return AAAA;
750 #else
751  throw Exceptions::RuntimeError("Epetra (Err. 2)");
752 #endif
753  break;
754  }
755  default:
756  throw Exceptions::RuntimeError("Only Epetra and Tpetra matrices can be transposed.");
757  break;
758  }
759 
760  return Teuchos::null;
761  }
762 
765  static RCP<Xpetra::MultiVector<double,LocalOrdinal,GlobalOrdinal,Node> > ExtractCoordinatesFromParameterList(ParameterList& paramList) {
766  RCP<Xpetra::MultiVector<double,LocalOrdinal,GlobalOrdinal,Node> > coordinates = Teuchos::null;
767 
768  // check whether coordinates are contained in parameter list
769  if(paramList.isParameter ("Coordinates") == false)
770  return coordinates;
771 
772  #if defined(HAVE_MUELU_TPETRA)
773  #if ( defined(EPETRA_HAVE_OMP) && defined(HAVE_TPETRA_INST_OPENMP) && defined(HAVE_TPETRA_INST_INT_INT)) || \
774  (!defined(EPETRA_HAVE_OMP) && defined(HAVE_TPETRA_INST_SERIAL) && defined(HAVE_TPETRA_INST_INT_INT))
775 
776  // define Tpetra::MultiVector type with Scalar=float only if
777  // * ETI is turned off, since then the compiler will instantiate it automatically OR
778  // * Tpetra is instantiated on Scalar=float
779  #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT)
781  RCP<tfMV> floatCoords = Teuchos::null;
782  #endif
783 
784  // define Tpetra::MultiVector type with Scalar=double only if
785  // * ETI is turned off, since then the compiler will instantiate it automatically OR
786  // * Tpetra is instantiated on Scalar=double
788  RCP<tdMV> doubleCoords = Teuchos::null;
789  if (paramList.isType<RCP<tdMV> >("Coordinates")) {
790  // Coordinates are stored as a double vector
791  doubleCoords = paramList.get<RCP<tdMV> >("Coordinates");
792  paramList.remove("Coordinates");
793  }
794  #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT)
795  else if (paramList.isType<RCP<tfMV> >("Coordinates")) {
796  // check if coordinates are stored as a float vector
797  floatCoords = paramList.get<RCP<tfMV> >("Coordinates");
798  paramList.remove("Coordinates");
799  doubleCoords = rcp(new tdMV(floatCoords->getMap(), floatCoords->getNumVectors()));
800  deep_copy(*doubleCoords, *floatCoords);
801  }
802  #endif
803  // We have the coordinates in a Tpetra double vector
804  if(doubleCoords != Teuchos::null) {
805  coordinates = Teuchos::rcp(new Xpetra::TpetraMultiVector<double, LocalOrdinal, GlobalOrdinal, Node>(doubleCoords));
806  TEUCHOS_TEST_FOR_EXCEPT(doubleCoords->getNumVectors() != coordinates->getNumVectors());
807  }
808  #endif // Tpetra instantiated on GO=int and EpetraNode
809  #endif // endif HAVE_TPETRA
810 
811  #if defined(HAVE_MUELU_EPETRA)
812  RCP<Epetra_MultiVector> doubleEpCoords;
813  if (paramList.isType<RCP<Epetra_MultiVector> >("Coordinates")) {
814  doubleEpCoords = paramList.get<RCP<Epetra_MultiVector> >("Coordinates");
815  paramList.remove("Coordinates");
816  RCP<Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> > epCoordinates = Teuchos::rcp(new Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node>(doubleEpCoords));
817  coordinates = rcp_dynamic_cast<Xpetra::MultiVector<double,LocalOrdinal,GlobalOrdinal,Node> >(epCoordinates);
818  TEUCHOS_TEST_FOR_EXCEPT(doubleEpCoords->NumVectors() != Teuchos::as<int>(coordinates->getNumVectors()));
819  }
820  #endif
821  TEUCHOS_TEST_FOR_EXCEPT(Teuchos::is_null(coordinates));
822  return coordinates;
823  }
824 
825  }; // class Utilities (specialization SC=double LO=GO=int)
826 
827 #endif // HAVE_MUELU_EPETRA
828 
829 
830 
835  long ExtractNonSerializableData(const Teuchos::ParameterList& inList, Teuchos::ParameterList& serialList, Teuchos::ParameterList& nonSerialList);
836 
837 
841  void TokenizeStringAndStripWhiteSpace(const std::string & stream, std::vector<std::string> & tokenList, const char* token = ",");
842 
845  bool IsParamMuemexVariable(const std::string& name);
846 
847 #ifdef HAVE_MUELU_EPETRA
848 
852  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
853  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
854  EpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Epetra_CrsMatrix>& A) {
858 
859  RCP<XCrsMatrix> Atmp = rcp(new XECrsMatrix(A));
860  return rcp(new XCrsMatrixWrap(Atmp));
861  }
862 
867  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
868  RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
869  EpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Epetra_MultiVector>& V) {
871  }
872 #endif
873 
874 #ifdef HAVE_MUELU_TPETRA
875 
879  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
880  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
885 
886  RCP<XCrsMatrix> Atmp = rcp(new XTCrsMatrix(Atpetra));
887  return rcp(new XCrsMatrixWrap(Atmp));
888  }
889 
894  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
895  RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
898  }
899 #endif
900 
902  template<class T>
903  std::string toString(const T& what) {
904  std::ostringstream buf;
905  buf << what;
906  return buf.str();
907  }
908 
909 #ifdef HAVE_MUELU_EPETRA
910 
914  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
915  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
916  EpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Epetra_CrsMatrix>& A);
917 
922  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
923  RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
924  EpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Epetra_MultiVector>& V);
925 #endif
926 
927 #ifdef HAVE_MUELU_TPETRA
928 
932  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
933  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
935 
940  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
941  RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
943 #endif
944 
945 } //namespace MueLu
946 
947 #define MUELU_UTILITIES_SHORT
948 #endif // MUELU_UTILITIES_DECL_HPP
static Teuchos::ArrayRCP< Scalar > GetMatrixDiagonal(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
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
static void SetRandomSeed(const Teuchos::Comm< int > &comm)
static Teuchos::ScalarTraits< Scalar >::magnitudeType Distance2(const MultiVector &v, LocalOrdinal i0, LocalOrdinal i1)
static RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > GetMatrixDiagonalInverse(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps() *100)
Exception indicating invalid cast attempted.
static Teuchos::RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > GetLumpedMatrixDiagonal(Teuchos::RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > A)
static void SetRandomSeed(const Teuchos::Comm< int > &comm)
static const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > Map2TpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
static const Epetra_MultiVector & MV2EpetraMV(const MultiVector &vec)
size_t getNodeMaxNumRowEntries() const
RCP< Epetra_CrsMatrix > getEpetra_CrsMatrixNonConst() 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.
RCP< Xpetra::Matrix< SC, LO, GO, NO > > EpetraCrs_To_XpetraMatrix(const Teuchos::RCP< Epetra_CrsMatrix > &A)
std::string toString(const T &what)
Little helper function to convert non-string types to strings.
RCP< Xpetra::MultiVector< SC, LO, GO, NO > > TpetraCrs_To_XpetraMultiVector(const Teuchos::RCP< Tpetra::MultiVector< SC, LO, GO, NO > > &Vtpetra)
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)
Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > Map
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())
virtual const RCP< const Map > & getRowMap() const
static RCP< Epetra_CrsMatrix > Op2NonConstEpetraCrs(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
void resumeFill(const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
static const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2TpetraCrs(const Matrix &Op)
static Scalar PowerMethod(const Matrix &A, bool scaleByDiag=true, LocalOrdinal niters=10, Magnitude tolerance=1e-2, bool verbose=false, unsigned int seed=123)
Simple transpose for Tpetra::CrsMatrix types.
static Scalar PowerMethod(const Matrix &A, bool scaleByDiag=true, LocalOrdinal niters=10, Magnitude tolerance=1e-2, bool verbose=false, unsigned int seed=123)
Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > CrsMatrix
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.
static Teuchos::ScalarTraits< Scalar >::magnitudeType Distance2(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &v, LocalOrdinal i0, LocalOrdinal i1)
static Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2NonConstTpetraCrs(Matrix &Op)
Teuchos::ScalarTraits< Scalar >::magnitudeType Magnitude
size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
Xpetra::CrsMatrixWrap< Scalar, LocalOrdinal, GlobalOrdinal, Node > CrsMatrixWrap
static RCP< const Epetra_CrsMatrix > Op2EpetraCrs(RCP< const Matrix > Op)
Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > Vector
void deep_copy(MultiVector< DS, DL, DG, DN, dstClassic > &dst, const MultiVector< SS, SL, SG, SN, srcClassic > &src)
static Teuchos::ArrayRCP< Scalar > GetLumpedMatrixDiagonal(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
static Teuchos::RCP< Vector > GetInverse(Teuchos::RCP< const Vector > v, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps() *100, Scalar tolReplacement=Teuchos::ScalarTraits< Scalar >::zero())
Return vector containing inverse of input vector.
Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > Matrix
static RCP< const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraRow(RCP< const Matrix > Op)
static Teuchos::ScalarTraits< Scalar >::magnitudeType Distance2(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &v, LocalOrdinal i0, LocalOrdinal i1)
Squared distance between two rows in a multivector.
static RCP< const Epetra_MultiVector > MV2EpetraMV(RCP< MultiVector > const vec)
Helper utility to pull out the underlying Epetra objects from an Xpetra object.
Namespace for MueLu classes and methods.
static const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > Map2TpetraMap(const Map &map)
Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > MultiVector
bool IsParamMuemexVariable(const std::string &name)
static RCP< Xpetra::MultiVector< double, LocalOrdinal, GlobalOrdinal, Node > > ExtractCoordinatesFromParameterList(ParameterList &paramList)
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec)
static Teuchos::RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > GetLumpedMatrixDiagonal(Teuchos::RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > A)
static void MyOldScaleMatrix_Tpetra(Matrix &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
static Scalar Frobenius(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B)
static RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraCrs(RCP< const Matrix > Op)
virtual void fillComplete(const RCP< const Map > &domainMap, const RCP< const Map > &rangeMap, const RCP< ParameterList > &params=null)=0
RCP< Epetra_MultiVector > getEpetra_MultiVector() const
static RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Residual(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &RHS)
static RCP< Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraRow(RCP< Matrix > Op)
static RCP< Epetra_MultiVector > MV2NonConstEpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec)
static RCP< Matrix > Transpose(Matrix &Op, bool optimizeTranspose=false, const std::string &label=std::string())
Transpose a Xpetra::Matrix.
static Epetra_MultiVector & MV2NonConstEpetraMV(MultiVector &vec)
bool isFillComplete() const
static RCP< Epetra_MultiVector > MV2NonConstEpetraMV(RCP< MultiVector > vec)
static RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraCrs(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static Teuchos::ArrayRCP< Scalar > GetMatrixDiagonal(const Matrix &A)
RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrixNonConst() const
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV(RCP< MultiVector > vec)
static Teuchos::ArrayRCP< Scalar > GetLumpedMatrixDiagonal(const Matrix &A)
static RCP< Xpetra::MultiVector< double, LocalOrdinal, GlobalOrdinal, Node > > ExtractCoordinatesFromParameterList(ParameterList &paramList)
Extract coordinates from parameter list and return them in a Xpetra::MultiVector. ...
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV2(Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
static const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & MV2TpetraMV(const MultiVector &vec)
RCP< CrsMatrix > getCrsMatrix() const
static Teuchos::ArrayRCP< const bool > DetectDirichletRows(const Matrix &A, const Magnitude &tol=Teuchos::ScalarTraits< Scalar >::zero())
static RCP< const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraRow(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
Xpetra::EpetraCrsMatrixT< GlobalOrdinal, Node > EpetraCrsMatrix
Teuchos::RCP< crs_matrix_type > createTranspose()
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV2(MultiVector &vec)
static void PauseForDebugger()
static Teuchos::ArrayRCP< Scalar > GetLumpedMatrixDiagonal(const Matrix &A)
Extract Matrix Diagonal of lumped matrix.
void getLocalRowView(LocalOrdinal LocalRow, Teuchos::ArrayView< const LocalOrdinal > &indices, Teuchos::ArrayView< const Scalar > &values) const
static Scalar Frobenius(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B)
Frobenius inner product of two matrices.
RCP< Xpetra::CrsMatrixWrap< SC, LO, GO, NO > > Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap(RCP< Epetra_CrsMatrix > &epAB)
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_CrsMatrix() const
RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > TpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &Vtpetra)
static RCP< Teuchos::FancyOStream > MakeFancy(std::ostream &os)
static Teuchos::Array< Magnitude > ResidualNorm(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const MultiVector &X, const MultiVector &RHS)
static RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Crs2Op(RCP< Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static void MyOldScaleMatrix_Epetra(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
RCP< const Epetra_CrsMatrix > getEpetra_CrsMatrix() const
static RCP< const Epetra_CrsMatrix > Op2EpetraCrs(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static Teuchos::RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > GetInverse(Teuchos::RCP< const Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > v, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps() *100, Scalar tolReplacement=Teuchos::ScalarTraits< Scalar >::zero())
static RCP< MultiVector > Residual(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const MultiVector &X, const MultiVector &RHS)
static RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > GetMatrixOverlappedDiagonal(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
static Teuchos::Array< Magnitude > ResidualNorm(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &RHS)
virtual bool isLocallyIndexed() const=0
static RCP< Vector > GetMatrixDiagonalInverse(const Matrix &A, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps() *100)
static RCP< Matrix > Crs2Op(RCP< CrsMatrix > Op)
Xpetra::CrsMatrixWrap< Scalar, LocalOrdinal, GlobalOrdinal, Node > CrsMatrixWrap
static RCP< Vector > GetMatrixDiagonalInverse(const Matrix &A, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps() *100)
Extract Matrix Diagonal.
RCP< Xpetra::Matrix< SC, LO, GO, NO > > TpetraCrs_To_XpetraMatrix(const Teuchos::RCP< Tpetra::CrsMatrix< SC, LO, GO, NO > > &Atpetra)
static RCP< Vector > GetMatrixOverlappedDiagonal(const Matrix &A)
static RCP< Teuchos::FancyOStream > MakeFancy(std::ostream &os)
static Teuchos::ArrayRCP< const bool > DetectDirichletRows(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Magnitude &tol=Teuchos::ScalarTraits< Scalar >::zero())
Detect Dirichlet rows.
Teuchos::RCP< const map_type > getDomainMap() const
static RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraCrs(RCP< Matrix > Op)
virtual Teuchos::RCP< const Map > getRangeMap() const=0
static Teuchos::ArrayRCP< const bool > DetectDirichletRows(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Magnitude &tol=Teuchos::ScalarTraits< Scalar >::magnitude(0.))
bool IsView(const viewLabel_t viewLabel) const
Exception throws to report errors in the internal logical of the program.
static Teuchos::Array< Magnitude > ResidualNorm(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &RHS)
RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_MultiVector() const
Teuchos::RCP< const map_type > getRangeMap() const
static Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & MV2NonConstTpetraMV(MultiVector &vec)
static const Epetra_Map & Map2EpetraMap(const Map &map)
static Scalar Frobenius(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B)
void TokenizeStringAndStripWhiteSpace(const std::string &stream, std::vector< std::string > &tokenList, const char *delimChars)
static void MyOldScaleMatrix(Matrix &Op, const Teuchos::ArrayRCP< const Scalar > &scalingVector, bool doInverse=true, bool doFillComplete=true, bool doOptimizeStorage=true)
static RCP< Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraRow(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< Vector > GetMatrixOverlappedDiagonal(const Matrix &A)
Extract Overlapped Matrix Diagonal.
static RCP< Vector > GetInverse(Teuchos::RCP< const Vector > v, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps() *100, Scalar tolReplacement=Teuchos::ScalarTraits< Scalar >::zero())
static RCP< Teuchos::FancyOStream > MakeFancy(std::ostream &os)
static void MyOldScaleMatrix_Tpetra(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
RCP< Xpetra::MultiVector< SC, LO, GO, NO > > EpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP< Epetra_MultiVector > &V)
static const Epetra_CrsMatrix & Op2EpetraCrs(const Matrix &Op)
static RCP< MultiVector > Residual(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &RHS)
static Scalar PowerMethod(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, bool scaleByDiag=true, LocalOrdinal niters=10, Magnitude tolerance=1e-2, bool verbose=false, unsigned int seed=123)
static RCP< Epetra_CrsMatrix > Op2NonConstEpetraCrs(RCP< Matrix > Op)
Teuchos::ScalarTraits< Scalar >::magnitudeType Magnitude
static void SetRandomSeed(const Teuchos::Comm< int > &comm)
Set seed for random number generator.
MueLu utility class.
virtual Teuchos::RCP< const Map > getDomainMap() const=0
static Teuchos::ArrayRCP< Scalar > GetMatrixDiagonal(const Matrix &A)
Extract Matrix Diagonal.
static void MyOldScaleMatrix_Epetra(Matrix &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
long ExtractNonSerializableData(const Teuchos::ParameterList &inList, Teuchos::ParameterList &serialList, Teuchos::ParameterList &nonSerialList)
Xpetra::EpetraMultiVectorT< GlobalOrdinal, Node > EpetraMultiVector