MueLu Version of the Day
Loading...
Searching...
No Matches
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 <string>
50
51#include "MueLu_ConfigDefs.hpp"
52
53#include <Teuchos_DefaultComm.hpp>
54#include <Teuchos_ScalarTraits.hpp>
55#include <Teuchos_ParameterList.hpp>
56
57#include <Xpetra_TpetraBlockCrsMatrix_fwd.hpp>
58#include <Xpetra_TpetraOperator.hpp>
59#include <Xpetra_CrsMatrix_fwd.hpp>
60#include <Xpetra_CrsMatrixWrap.hpp>
61#include <Xpetra_Map_fwd.hpp>
62#include <Xpetra_Matrix_fwd.hpp>
63#include <Xpetra_MultiVector_fwd.hpp>
64#include <Xpetra_MultiVectorFactory_fwd.hpp>
65#include <Xpetra_Operator_fwd.hpp>
66#include <Xpetra_Vector_fwd.hpp>
67#include <Xpetra_VectorFactory_fwd.hpp>
68
69#include <Xpetra_MatrixMatrix.hpp>
70
71#ifdef HAVE_MUELU_EPETRA
72#include <Xpetra_EpetraCrsMatrix.hpp>
73
74// needed because of inlined function
75//TODO: remove inline function?
76#include <Xpetra_EpetraCrsMatrix_fwd.hpp>
77#include <Xpetra_CrsMatrixWrap_fwd.hpp>
78
79#endif
80
81#include "MueLu_Exceptions.hpp"
82
83#ifdef HAVE_MUELU_EPETRAEXT
86class Epetra_Vector;
88#endif
89
90#include <Tpetra_CrsMatrix.hpp>
91#include <Tpetra_BlockCrsMatrix.hpp>
92#include <Tpetra_BlockCrsMatrix_Helpers.hpp>
93#include <Tpetra_FECrsMatrix.hpp>
94#include <Tpetra_RowMatrixTransposer.hpp>
95#include <Tpetra_Map.hpp>
96#include <Tpetra_MultiVector.hpp>
97#include <Tpetra_FEMultiVector.hpp>
98#include <Xpetra_TpetraCrsMatrix_fwd.hpp>
99#include <Xpetra_TpetraMultiVector_fwd.hpp>
100
101#include <MueLu_UtilitiesBase.hpp>
102
103
104namespace MueLu {
105
106#ifdef HAVE_MUELU_EPETRA
107 //defined after Utilities class
108 template<typename SC,typename LO,typename GO,typename NO>
109 RCP<Xpetra::CrsMatrixWrap<SC,LO,GO,NO> >
110 Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap(RCP<Epetra_CrsMatrix> &epAB);
111
112 template<typename SC,typename LO,typename GO,typename NO>
113 RCP<Xpetra::Matrix<SC, LO, GO, NO> >
114 EpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Epetra_CrsMatrix>& A);
115
116 template<typename SC,typename LO,typename GO,typename NO>
117 RCP<Xpetra::MultiVector<SC, LO, GO, NO> >
118 EpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Epetra_MultiVector>& V);
119#endif
120
121 template<typename SC,typename LO,typename GO,typename NO>
122 RCP<Xpetra::Matrix<SC, LO, GO, NO> >
123 TpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Tpetra::CrsMatrix<SC, LO, GO, NO> >& Atpetra);
124
125 template<typename SC,typename LO,typename GO,typename NO>
126 RCP<Xpetra::Matrix<SC, LO, GO, NO> >
127 TpetraFECrs_To_XpetraMatrix(const Teuchos::RCP<Tpetra::FECrsMatrix<SC, LO, GO, NO> >& Atpetra);
128
129 template<typename SC,typename LO,typename GO,typename NO>
130 RCP<Xpetra::MultiVector<SC, LO, GO, NO> >
131 TpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Tpetra::MultiVector<SC, LO, GO, NO> >& Vtpetra);
132
133 template<typename SC,typename LO,typename GO,typename NO>
134 RCP<Xpetra::MultiVector<SC, LO, GO, NO> >
135 TpetraFEMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Tpetra::FEMultiVector<SC, LO, GO, NO> >& Vtpetra);
136
137 template<typename SC,typename LO,typename GO,typename NO>
138 void leftRghtDofScalingWithinNode(const Xpetra::Matrix<SC,LO,GO,NO> & Atpetra, size_t blkSize, size_t nSweeps, Teuchos::ArrayRCP<SC> & rowScaling, Teuchos::ArrayRCP<SC> & colScaling);
139
147 template <class Scalar,
150 class Node = DefaultNode>
151 class Utilities : public UtilitiesBase<Scalar, LocalOrdinal, GlobalOrdinal, Node> {
152#undef MUELU_UTILITIES_SHORT
153 #include "MueLu_UseShortNames.hpp"
154
155 public:
156 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType Magnitude;
157
158#ifdef HAVE_MUELU_EPETRA
160 // @{
161 static RCP<const Epetra_MultiVector> MV2EpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > const vec);
162 static RCP< Epetra_MultiVector> MV2NonConstEpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > vec);
163
164 static const Epetra_MultiVector& MV2EpetraMV(const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
165 static Epetra_MultiVector& MV2NonConstEpetraMV(Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
166
167 static RCP<const Epetra_CrsMatrix> Op2EpetraCrs(RCP<const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
168 static RCP< Epetra_CrsMatrix> Op2NonConstEpetraCrs(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
169
170 static const Epetra_CrsMatrix& Op2EpetraCrs(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
171 static Epetra_CrsMatrix& Op2NonConstEpetraCrs(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
172
173 static const Epetra_Map& Map2EpetraMap(const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node>& map);
174 // @}
175#endif
176
178 static RCP<const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2TpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > const vec);
179 static RCP< Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2NonConstTpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > vec);
180 static RCP< Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2NonConstTpetraMV2(Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
181
182 static const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& MV2TpetraMV(const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
183 static Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& MV2NonConstTpetraMV(Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
184
185 static RCP<const Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2TpetraCrs(RCP<const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
186 static RCP< Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2NonConstTpetraCrs(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
187
188 static const Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2TpetraCrs(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
189 static Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2NonConstTpetraCrs(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
190
191 static RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2TpetraBlockCrs(RCP<const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
192 static RCP< Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2NonConstTpetraBlockCrs(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
193
194 static const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2TpetraBlockCrs(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
195 static Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2NonConstTpetraBlockCrs(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
196
197
198
199 static RCP<const Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2TpetraRow(RCP<const Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
200 static RCP< Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2NonConstTpetraRow(RCP<Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
201
202
203 static const RCP<const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > Map2TpetraMap(const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node>& map);
204
205 static void MyOldScaleMatrix(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Teuchos::ArrayRCP<const Scalar>& scalingVector, bool doInverse = true,
206 bool doFillComplete = true, bool doOptimizeStorage = true);
207
208 static void MyOldScaleMatrix_Epetra(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Teuchos::ArrayRCP<Scalar>& scalingVector,
209 bool doFillComplete, bool doOptimizeStorage);
210 static void MyOldScaleMatrix_Tpetra(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Teuchos::ArrayRCP<Scalar>& scalingVector,
211 bool doFillComplete, bool doOptimizeStorage);
212
213 static RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Transpose(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, bool optimizeTranspose = false,const std::string & label = std::string(),const Teuchos::RCP<Teuchos::ParameterList> &params=Teuchos::null);
214
215 static RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > RealValuedToScalarMultiVector(RCP<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::coordinateType,LocalOrdinal,GlobalOrdinal,Node> > X);
216
217 static RCP<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LocalOrdinal,GlobalOrdinal,Node> > ExtractCoordinatesFromParameterList(ParameterList& paramList);
218
219 }; // class Utilities
220
222
223#ifdef HAVE_MUELU_EPETRA
233 template <>
234 class Utilities<double,int,int,Xpetra::EpetraNode> : public UtilitiesBase<double,int,int,Xpetra::EpetraNode> {
235 public:
236 typedef double Scalar;
237 typedef int LocalOrdinal;
238 typedef int GlobalOrdinal;
239 typedef Xpetra::EpetraNode Node;
240 typedef Teuchos::ScalarTraits<Scalar>::magnitudeType Magnitude;
241#undef MUELU_UTILITIES_SHORT
243
244 private:
245 using EpetraMap = Xpetra::EpetraMapT<GlobalOrdinal,Node>;
246 using EpetraMultiVector = Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node>;
247 // using EpetraCrsMatrix = Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>;
248 public:
249
251 // @{
252 static RCP<const Epetra_MultiVector> MV2EpetraMV(RCP<MultiVector> const vec) {
253 RCP<const EpetraMultiVector > tmpVec = rcp_dynamic_cast<EpetraMultiVector>(vec);
254 if (tmpVec == Teuchos::null)
255 throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
256 return tmpVec->getEpetra_MultiVector();
257 }
258 static RCP< Epetra_MultiVector> MV2NonConstEpetraMV(RCP<MultiVector> vec) {
259 RCP<const EpetraMultiVector> tmpVec = rcp_dynamic_cast<EpetraMultiVector>(vec);
260 if (tmpVec == Teuchos::null)
261 throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
262 return tmpVec->getEpetra_MultiVector();
263 }
264
265 static const Epetra_MultiVector& MV2EpetraMV(const MultiVector& vec) {
266 const EpetraMultiVector& tmpVec = dynamic_cast<const EpetraMultiVector&>(vec);
267 return *(tmpVec.getEpetra_MultiVector());
268 }
269 static Epetra_MultiVector& MV2NonConstEpetraMV(MultiVector& vec) {
270 const EpetraMultiVector& tmpVec = dynamic_cast<const EpetraMultiVector&>(vec);
271 return *(tmpVec.getEpetra_MultiVector());
272 }
273
274 static RCP<const Epetra_CrsMatrix> Op2EpetraCrs(RCP<const Matrix> Op) {
275 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
276 if (crsOp == Teuchos::null)
277 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
278 const RCP<const EpetraCrsMatrix>& tmp_ECrsMtx = rcp_dynamic_cast<const EpetraCrsMatrix>(crsOp->getCrsMatrix());
279 if (tmp_ECrsMtx == Teuchos::null)
280 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
281 return tmp_ECrsMtx->getEpetra_CrsMatrix();
282 }
283 static RCP< Epetra_CrsMatrix> Op2NonConstEpetraCrs(RCP<Matrix> Op) {
284 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
285 if (crsOp == Teuchos::null)
286 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
287 const RCP<const EpetraCrsMatrix> &tmp_ECrsMtx = rcp_dynamic_cast<const EpetraCrsMatrix>(crsOp->getCrsMatrix());
288 if (tmp_ECrsMtx == Teuchos::null)
289 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
290 return tmp_ECrsMtx->getEpetra_CrsMatrixNonConst();
291 }
292
293 static const Epetra_CrsMatrix& Op2EpetraCrs(const Matrix& Op) {
294 try {
295 const CrsMatrixWrap& crsOp = dynamic_cast<const CrsMatrixWrap&>(Op);
296 try {
297 const EpetraCrsMatrix& tmp_ECrsMtx = dynamic_cast<const EpetraCrsMatrix&>(*crsOp.getCrsMatrix());
298 return *tmp_ECrsMtx.getEpetra_CrsMatrix();
299 } catch (std::bad_cast&) {
300 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
301 }
302 } catch (std::bad_cast&) {
303 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
304 }
305 }
307 try {
308 CrsMatrixWrap& crsOp = dynamic_cast<CrsMatrixWrap&>(Op);
309 try {
310 EpetraCrsMatrix& tmp_ECrsMtx = dynamic_cast<EpetraCrsMatrix&>(*crsOp.getCrsMatrix());
311 return *tmp_ECrsMtx.getEpetra_CrsMatrixNonConst();
312 } catch (std::bad_cast&) {
313 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
314 }
315 } catch (std::bad_cast&) {
316 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
317 }
318 }
319
320 static const Epetra_Map& Map2EpetraMap(const Map& map) {
321 RCP<const EpetraMap> xeMap = rcp_dynamic_cast<const EpetraMap>(rcpFromRef(map));
322 if (xeMap == Teuchos::null)
323 throw Exceptions::BadCast("Utilities::Map2EpetraMap : Cast from Xpetra::Map to Xpetra::EpetraMap failed");
324 return xeMap->getEpetra_Map();
325 }
326 // @}
327
329 static RCP<const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2TpetraMV(RCP<MultiVector> const vec) {
330#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
331 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
332 throw Exceptions::RuntimeError("MV2TpetraMV: Tpetra has not been compiled with support for LO=GO=int.");
333#else
334 RCP<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(vec);
335 if (tmpVec == Teuchos::null)
336 throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
337 return tmpVec->getTpetra_MultiVector();
338#endif
339 }
340 static RCP< Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2NonConstTpetraMV(RCP<MultiVector> vec) {
341#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
342 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
343 throw Exceptions::RuntimeError("MV2NonConstTpetraMV: Tpetra has not been compiled with support for LO=GO=int.");
344#else
345 RCP<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(vec);
346 if (tmpVec == Teuchos::null)
347 throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
348 return tmpVec->getTpetra_MultiVector();
349#endif
350
351 }
352 static RCP< Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2NonConstTpetraMV2(MultiVector& vec) {
353#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
354 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
355 throw Exceptions::RuntimeError("MV2NonConstTpetraMV2: Tpetra has not been compiled with support for LO=GO=int.");
356#else
357 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec = dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(vec);
358 return tmpVec.getTpetra_MultiVector();
359#endif
360 }
361
362 static const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& MV2TpetraMV(const MultiVector& vec) {
363#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
364 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
365 throw Exceptions::RuntimeError("MV2TpetraMV: Tpetra has not been compiled with support for LO=GO=int.");
366#else
367 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec = dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(vec);
368 return *(tmpVec.getTpetra_MultiVector());
369#endif
370 }
371 static Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& MV2NonConstTpetraMV(MultiVector& vec) {
372#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
373 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
374 throw Exceptions::RuntimeError("MV2NonConstTpetraMV: Tpetra has not been compiled with support for LO=GO=int.");
375#else
376 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec = dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(vec);
377 return *(tmpVec.getTpetra_MultiVector());
378#endif
379 }
380
381 static RCP<const Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2TpetraCrs(RCP<const Matrix> Op) {
382#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
383 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
384 throw Exceptions::RuntimeError("Op2TpetraCrs: Tpetra has not been compiled with support for LO=GO=int.");
385#else
386 // Get the underlying Tpetra Mtx
387 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
388 if (crsOp == Teuchos::null)
389 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
390 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
391 if (tmp_ECrsMtx == Teuchos::null)
392 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
393 return tmp_ECrsMtx->getTpetra_CrsMatrix();
394#endif
395 }
396 static RCP< Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2NonConstTpetraCrs(RCP<Matrix> Op){
397#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
398 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
399 throw Exceptions::RuntimeError("Op2NonConstTpetraCrs: Tpetra has not been compiled with support for LO=GO=int.");
400#else
401 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
402 if (crsOp == Teuchos::null)
403 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
404 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
405 if (tmp_ECrsMtx == Teuchos::null)
406 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
407 return tmp_ECrsMtx->getTpetra_CrsMatrixNonConst();
408#endif
409 };
410
411 static const Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2TpetraCrs(const Matrix& Op) {
412#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
413 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
414 throw Exceptions::RuntimeError("Op2TpetraCrs: Tpetra has not been compiled with support for LO=GO=int.");
415#else
416 try {
417 const CrsMatrixWrap& crsOp = dynamic_cast<const CrsMatrixWrap&>(Op);
418 try {
419 const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx = dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(*crsOp.getCrsMatrix());
420 return *tmp_ECrsMtx.getTpetra_CrsMatrix();
421 } catch (std::bad_cast&) {
422 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
423 }
424 } catch (std::bad_cast&) {
425 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
426 }
427#endif
428 }
429 static Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2NonConstTpetraCrs(Matrix& Op) {
430#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
431 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
432 throw Exceptions::RuntimeError("Op2NonConstTpetraCrs: Tpetra has not been compiled with support for LO=GO=int.");
433#else
434 try {
435 CrsMatrixWrap& crsOp = dynamic_cast<CrsMatrixWrap&>(Op);
436 try {
437 Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx = dynamic_cast<Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(*crsOp.getCrsMatrix());
438 return *tmp_ECrsMtx.getTpetra_CrsMatrixNonConst();
439 } catch (std::bad_cast&) {
440 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
441 }
442 } catch (std::bad_cast&) {
443 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
444 }
445#endif
446 }
447
448
449 static RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2TpetraBlockCrs(RCP<const Matrix> Op) {
450#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
451 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
452 throw Exceptions::RuntimeError("Op2TpetraBlockCrs: Tpetra has not been compiled with support for LO=GO=int.");
453#else
454 // Get the underlying Tpetra Mtx
455 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
456 if (crsOp == Teuchos::null)
457 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
458 const RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
459 if (tmp_ECrsMtx == Teuchos::null)
460 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
461 return tmp_ECrsMtx->getTpetra_BlockCrsMatrix();
462#endif
463 }
464
465 static RCP< Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2NonConstTpetraBlockCrs(RCP<Matrix> Op){
466#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
467 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
468 throw Exceptions::RuntimeError("Op2NonConstTpetraBlockCrs: Tpetra has not been compiled with support for LO=GO=int.");
469#else
470 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(Op);
471 if (crsOp == Teuchos::null)
472 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
473 const RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
474 if (tmp_ECrsMtx == Teuchos::null)
475 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
476 return tmp_ECrsMtx->getTpetra_BlockCrsMatrixNonConst();
477#endif
478 };
479
480 static const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2TpetraBlockCrs(const Matrix& Op) {
481#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
482 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
483 throw Exceptions::RuntimeError("Op2TpetraBlockCrs: Tpetra has not been compiled with support for LO=GO=int.");
484#else
485 try {
486 const CrsMatrixWrap& crsOp = dynamic_cast<const CrsMatrixWrap&>(Op);
487 try {
488 const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx = dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(*crsOp.getCrsMatrix());
489 return *tmp_ECrsMtx.getTpetra_BlockCrsMatrix();
490 } catch (std::bad_cast&) {
491 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
492 }
493 } catch (std::bad_cast&) {
494 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
495 }
496#endif
497 }
498 static Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2NonConstTpetraBlockCrs(Matrix& Op) {
499#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
500 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
501 throw Exceptions::RuntimeError("Op2NonConstTpetraCrs: Tpetra has not been compiled with support for LO=GO=int.");
502#else
503 try {
504 CrsMatrixWrap& crsOp = dynamic_cast<CrsMatrixWrap&>(Op);
505 try {
506 Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx = dynamic_cast<Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(*crsOp.getCrsMatrix());
507 return *tmp_ECrsMtx.getTpetra_BlockCrsMatrixNonConst();
508 } catch (std::bad_cast&) {
509 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
510 }
511 } catch (std::bad_cast&) {
512 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
513 }
514#endif
515 }
516
517
518 static RCP<const Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2TpetraRow(RCP<const Operator> Op) {
519#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
520 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
521 throw Exceptions::RuntimeError("Op2TpetraRow: Tpetra has not been compiled with support for LO=GO=int.");
522#else
523 RCP<const Matrix> mat = rcp_dynamic_cast<const Matrix>(Op);
524 if (!mat.is_null()) {
525 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(mat);
526 if (crsOp == Teuchos::null)
527 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
528
529 RCP<const CrsMatrix> crsMat = crsOp->getCrsMatrix();
530 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_Crs = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
531 RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_BlockCrs;
532 if(!tmp_Crs.is_null()) {
533 return tmp_Crs->getTpetra_CrsMatrixNonConst();
534 }
535 else {
536 tmp_BlockCrs= rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
537 if (tmp_BlockCrs.is_null())
538 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
539 return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
540 }
541 } else {
542 RCP<const TpetraOperator> tpOp = rcp_dynamic_cast<const TpetraOperator>(Op, true);
543 RCP<const Tpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node>> tOp = tpOp->getOperatorConst();
544 RCP<const Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tRow = rcp_dynamic_cast<const Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(tOp, true);
545 return tRow;
546 }
547#endif
548 }
549
550 static RCP< Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2NonConstTpetraRow(RCP<Operator> Op) {
551#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
552 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
553 throw Exceptions::RuntimeError("Op2NonConstTpetraRow: Tpetra has not been compiled with support for LO=GO=int.");
554#else
555 RCP<Matrix> mat = rcp_dynamic_cast<Matrix>(Op);
556 if (!mat.is_null()) {
557 RCP<const CrsMatrixWrap> crsOp = rcp_dynamic_cast<const CrsMatrixWrap>(mat);
558 if (crsOp == Teuchos::null)
559 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
560
561 RCP<const CrsMatrix> crsMat = crsOp->getCrsMatrix();
562 const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_Crs = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
563 RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_BlockCrs;
564 if(!tmp_Crs.is_null()) {
565 return tmp_Crs->getTpetra_CrsMatrixNonConst();
566 }
567 else {
568 tmp_BlockCrs= rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
569 if (tmp_BlockCrs.is_null())
570 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
571 return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
572 }
573 } else {
574 RCP<TpetraOperator> tpOp = rcp_dynamic_cast<TpetraOperator>(Op, true);
575 RCP<Tpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node>> tOp = tpOp->getOperator();
576 RCP<Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tRow = rcp_dynamic_cast<Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(tOp, true);
577 return tRow;
578 }
579#endif
580 };
581
582
583 static const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > Map2TpetraMap(const Map& map) {
584#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
585 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
586 throw Exceptions::RuntimeError("Map2TpetraMap: Tpetra has not been compiled with support for LO=GO=int.");
587#else
588 const RCP<const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node>>& tmp_TMap = rcp_dynamic_cast<const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node> >(rcpFromRef(map));
589 if (tmp_TMap == Teuchos::null)
590 throw Exceptions::BadCast("Utilities::Map2TpetraMap : Cast from Xpetra::Map to Xpetra::TpetraMap failed");
591 return tmp_TMap->getTpetra_Map();
592#endif
593 };
594
595
596 static void MyOldScaleMatrix(Matrix& Op, const Teuchos::ArrayRCP<const Scalar>& scalingVector, bool doInverse = true,
597 bool doFillComplete = true, bool doOptimizeStorage = true) {
598 Scalar one = Teuchos::ScalarTraits<Scalar>::one();
599 Teuchos::ArrayRCP<Scalar> sv(scalingVector.size());
600 if (doInverse) {
601 for (int i = 0; i < scalingVector.size(); ++i)
602 sv[i] = one / scalingVector[i];
603 } else {
604 for (int i = 0; i < scalingVector.size(); ++i)
605 sv[i] = scalingVector[i];
606 }
607
608 switch (Op.getRowMap()->lib()) {
609 case Xpetra::UseTpetra:
610 MyOldScaleMatrix_Tpetra(Op, sv, doFillComplete, doOptimizeStorage);
611 break;
612
613 case Xpetra::UseEpetra:
614 MyOldScaleMatrix_Epetra(Op, sv, doFillComplete, doOptimizeStorage);
615 break;
616
617 default:
618 throw Exceptions::RuntimeError("Only Epetra and Tpetra matrices can be scaled.");
619 }
620 }
621
622 // TODO This is the <double,int,int> specialization
623 static void MyOldScaleMatrix_Tpetra(Matrix& Op, const Teuchos::ArrayRCP<Scalar>& scalingVector,
624 bool doFillComplete, bool doOptimizeStorage) {
625#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
626 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
627 throw Exceptions::RuntimeError("Matrix scaling is not possible because Tpetra has not been compiled with support for LO=GO=int.");
628#else
629 try {
630 Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& tpOp = Op2NonConstTpetraCrs(Op);
631
632 const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rowMap = tpOp.getRowMap();
633 const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > domainMap = tpOp.getDomainMap();
634 const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rangeMap = tpOp.getRangeMap();
635
636 size_t maxRowSize = tpOp.getLocalMaxNumRowEntries();
637 if (maxRowSize == Teuchos::as<size_t>(-1)) // hasn't been determined yet
638 maxRowSize = 20;
639
640 std::vector<Scalar> scaledVals(maxRowSize);
641 if (tpOp.isFillComplete())
642 tpOp.resumeFill();
643
644 if (Op.isLocallyIndexed() == true) {
645 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::local_inds_host_view_type cols;
646 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::values_host_view_type vals;
647 for (size_t i = 0; i < rowMap->getLocalNumElements(); ++i) {
648 tpOp.getLocalRowView(i, cols, vals);
649 size_t nnz = tpOp.getNumEntriesInLocalRow(i);
650 if (nnz > maxRowSize) {
651 maxRowSize = nnz;
652 scaledVals.resize(maxRowSize);
653 }
654 for (size_t j = 0; j < nnz; ++j)
655 scaledVals[j] = vals[j]*scalingVector[i];
656
657 if (nnz > 0) {
658 Teuchos::ArrayView<const LocalOrdinal> cols_view(cols.data(), nnz);
659 Teuchos::ArrayView<const Scalar> valview(&scaledVals[0], nnz);
660 tpOp.replaceLocalValues(i, cols_view, valview);
661 }
662 } //for (size_t i=0; ...
663
664 } else {
665 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::local_inds_host_view_type cols;
666 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::values_host_view_type vals;
667
668 for (size_t i = 0; i < rowMap->getLocalNumElements(); ++i) {
669 GlobalOrdinal gid = rowMap->getGlobalElement(i);
670 tpOp.getGlobalRowView(gid, cols, vals);
671 size_t nnz = tpOp.getNumEntriesInGlobalRow(gid);
672 if (nnz > maxRowSize) {
673 maxRowSize = nnz;
674 scaledVals.resize(maxRowSize);
675 }
676 // FIXME FIXME FIXME FIXME FIXME FIXME
677 for (size_t j = 0; j < nnz; ++j)
678 scaledVals[j] = vals[j]*scalingVector[i]; //FIXME i or gid?
679
680 if (nnz > 0) {
681 Teuchos::ArrayView<const LocalOrdinal> cols_view(cols.data(), nnz);
682 Teuchos::ArrayView<const Scalar> valview(&scaledVals[0], nnz);
683 tpOp.replaceGlobalValues(gid, cols_view, valview);
684 }
685 } //for (size_t i=0; ...
686 }
687
688 if (doFillComplete) {
689 if (domainMap == Teuchos::null || rangeMap == Teuchos::null)
690 throw Exceptions::RuntimeError("In Utilities::Scaling: cannot fillComplete because the domain and/or range map hasn't been defined");
691
692 RCP<Teuchos::ParameterList> params = rcp(new Teuchos::ParameterList());
693 params->set("Optimize Storage", doOptimizeStorage);
694 params->set("No Nonlocal Changes", true);
695 Op.fillComplete(Op.getDomainMap(), Op.getRangeMap(), params);
696 }
697 } catch(...) {
698 throw Exceptions::RuntimeError("Only Tpetra::CrsMatrix types can be scaled (Err.1)");
699 }
700#endif
701 }
702
703 static void MyOldScaleMatrix_Epetra (Matrix& Op, const Teuchos::ArrayRCP<Scalar>& scalingVector, bool /* doFillComplete */, bool /* doOptimizeStorage */) {
704#ifdef HAVE_MUELU_EPETRA
705 try {
706 //const Epetra_CrsMatrix& epOp = Utilities<double,int,int>::Op2NonConstEpetraCrs(Op);
707 const Epetra_CrsMatrix& epOp = Op2NonConstEpetraCrs(Op);
708
709 Epetra_Map const &rowMap = epOp.RowMap();
710 int nnz;
711 double *vals;
712 int *cols;
713
714 for (int i = 0; i < rowMap.NumMyElements(); ++i) {
715 epOp.ExtractMyRowView(i, nnz, vals, cols);
716 for (int j = 0; j < nnz; ++j)
717 vals[j] *= scalingVector[i];
718 }
719
720 } catch (...){
721 throw Exceptions::RuntimeError("Only Epetra_CrsMatrix types can be scaled");
722 }
723#else
724 throw Exceptions::RuntimeError("Matrix scaling is not possible because Epetra has not been enabled.");
725#endif // HAVE_MUELU_EPETRA
726 }
727
733 static RCP<Matrix> Transpose(Matrix& Op, bool /* optimizeTranspose */ = false,const std::string & label = std::string(),const Teuchos::RCP<Teuchos::ParameterList> &params=Teuchos::null) {
734 switch (Op.getRowMap()->lib()) {
735 case Xpetra::UseTpetra: {
736#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
737 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
738 throw Exceptions::RuntimeError("Utilities::Transpose: Tpetra is not compiled with LO=GO=int. Add TPETRA_INST_INT_INT:BOOL=ON to your configuration!");
739#else
740 using Helpers = Xpetra::Helpers<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
741 /***************************************************************/
742 if(Helpers::isTpetraCrs(Op)) {
743 const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& tpetraOp = Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Op2TpetraCrs(Op);
744
745 // Compute the transpose A of the Tpetra matrix tpetraOp.
746 RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > A;
747 Tpetra::RowMatrixTransposer<Scalar, LocalOrdinal, GlobalOrdinal, Node> transposer(rcpFromRef(tpetraOp),label);
748
749 {
750 using Teuchos::ParameterList;
751 using Teuchos::rcp;
752 RCP<ParameterList> transposeParams = params.is_null () ?
753 rcp (new ParameterList) :
754 rcp (new ParameterList (*params));
755 transposeParams->set ("sort", false);
756 A = transposer.createTranspose(transposeParams);
757 }
758
759 RCP<Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AA = rcp(new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(A));
760 RCP<CrsMatrix> AAA = rcp_implicit_cast<CrsMatrix>(AA);
761 RCP<Matrix> AAAA = rcp( new CrsMatrixWrap(AAA));
762
763 if (Op.IsView("stridedMaps"))
764 AAAA->CreateView("stridedMaps", Teuchos::rcpFromRef(Op), true/*doTranspose*/);
765
766 return AAAA;
767 }
768 /***************************************************************/
769 else if(Helpers::isTpetraBlockCrs(Op)) {
770 using BCRS = Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
771 // using CRS = Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
773 RCP<BCRS> At;
774 {
775 Tpetra::BlockCrsMatrixTransposer<Scalar, LocalOrdinal, GlobalOrdinal, Node> transposer(rcpFromRef(tpetraOp),label);
776
777 using Teuchos::ParameterList;
778 using Teuchos::rcp;
779 RCP<ParameterList> transposeParams = params.is_null () ?
780 rcp (new ParameterList) :
781 rcp (new ParameterList (*params));
782 transposeParams->set ("sort", false);
783 At = transposer.createTranspose(transposeParams);
784 }
785
786 RCP<Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AA = rcp(new Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(At));
787 RCP<CrsMatrix> AAA = rcp_implicit_cast<CrsMatrix>(AA);
788 RCP<Matrix> AAAA = rcp( new CrsMatrixWrap(AAA));
789
790 if (Op.IsView("stridedMaps"))
791 AAAA->CreateView("stridedMaps", Teuchos::rcpFromRef(Op), true/*doTranspose*/);
792
793 return AAAA;
794
795 }
796 /***************************************************************/
797 else {
798 throw Exceptions::RuntimeError("Utilities::Transpose failed, perhaps because matrix is not a Crs or BlockCrs matrix");
799 }
800#endif
801 }
802 case Xpetra::UseEpetra:
803 {
804#if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT)
805 Teuchos::TimeMonitor tm(*Teuchos::TimeMonitor::getNewTimer("ZZ Entire Transpose"));
806 // Epetra case
809 Epetra_CrsMatrix * A = dynamic_cast<Epetra_CrsMatrix*>(&transposer(epetraOp));
810 transposer.ReleaseTranspose(); // So we can keep A in Muelu...
811
812 RCP<Epetra_CrsMatrix> rcpA(A);
813 RCP<EpetraCrsMatrix> AA = rcp(new EpetraCrsMatrix(rcpA));
814 RCP<CrsMatrix> AAA = rcp_implicit_cast<CrsMatrix>(AA);
815 RCP<Matrix> AAAA = rcp( new CrsMatrixWrap(AAA));
816
817 if (Op.IsView("stridedMaps"))
818 AAAA->CreateView("stridedMaps", Teuchos::rcpFromRef(Op), true/*doTranspose*/);
819
820 return AAAA;
821#else
822 throw Exceptions::RuntimeError("Epetra (Err. 2)");
823#endif
824 }
825 default:
826 throw Exceptions::RuntimeError("Only Epetra and Tpetra matrices can be transposed.");
827 }
828
829 TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
830 }
831
832 static RCP<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LocalOrdinal,GlobalOrdinal,Node> >
833 RealValuedToScalarMultiVector(RCP<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::coordinateType,LocalOrdinal,GlobalOrdinal,Node> > X) {
834 RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Xscalar = rcp_dynamic_cast<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(X, true);
835 return Xscalar;
836 }
837
840 static RCP<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LocalOrdinal,GlobalOrdinal,Node> > ExtractCoordinatesFromParameterList(ParameterList& paramList) {
841 RCP<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LocalOrdinal,GlobalOrdinal,Node> > coordinates = Teuchos::null;
842
843 // check whether coordinates are contained in parameter list
844 if(paramList.isParameter ("Coordinates") == false)
845 return coordinates;
846
847 #if ( defined(EPETRA_HAVE_OMP) && defined(HAVE_TPETRA_INST_OPENMP) && defined(HAVE_TPETRA_INST_INT_INT)) || \
848 (!defined(EPETRA_HAVE_OMP) && defined(HAVE_TPETRA_INST_SERIAL) && defined(HAVE_TPETRA_INST_INT_INT))
849
850 // define Tpetra::MultiVector type with Scalar=float only if
851 // * ETI is turned off, since then the compiler will instantiate it automatically OR
852 // * Tpetra is instantiated on Scalar=float
853 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT)
854 typedef Tpetra::MultiVector<float, LocalOrdinal, GlobalOrdinal, Node> tfMV;
855 RCP<tfMV> floatCoords = Teuchos::null;
856 #endif
857
858 // define Tpetra::MultiVector type with Scalar=double only if
859 // * ETI is turned off, since then the compiler will instantiate it automatically OR
860 // * Tpetra is instantiated on Scalar=double
861 typedef Tpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType, LocalOrdinal, GlobalOrdinal, Node> tdMV;
862 RCP<tdMV> doubleCoords = Teuchos::null;
863 if (paramList.isType<RCP<tdMV> >("Coordinates")) {
864 // Coordinates are stored as a double vector
865 doubleCoords = paramList.get<RCP<tdMV> >("Coordinates");
866 paramList.remove("Coordinates");
867 }
868 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT)
869 else if (paramList.isType<RCP<tfMV> >("Coordinates")) {
870 // check if coordinates are stored as a float vector
871 floatCoords = paramList.get<RCP<tfMV> >("Coordinates");
872 paramList.remove("Coordinates");
873 doubleCoords = rcp(new tdMV(floatCoords->getMap(), floatCoords->getNumVectors()));
874 deep_copy(*doubleCoords, *floatCoords);
875 }
876 #endif
877 // We have the coordinates in a Tpetra double vector
878 if(doubleCoords != Teuchos::null) {
879 coordinates = Teuchos::rcp(new Xpetra::TpetraMultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType, LocalOrdinal, GlobalOrdinal, Node>(doubleCoords));
880 TEUCHOS_TEST_FOR_EXCEPT(doubleCoords->getNumVectors() != coordinates->getNumVectors());
881 }
882 #endif // Tpetra instantiated on GO=int and EpetraNode
883
884 #if defined(HAVE_MUELU_EPETRA)
885 RCP<Epetra_MultiVector> doubleEpCoords;
886 if (paramList.isType<RCP<Epetra_MultiVector> >("Coordinates")) {
887 doubleEpCoords = paramList.get<RCP<Epetra_MultiVector> >("Coordinates");
888 paramList.remove("Coordinates");
889 RCP<Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> > epCoordinates = Teuchos::rcp(new Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node>(doubleEpCoords));
890 coordinates = rcp_dynamic_cast<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LocalOrdinal,GlobalOrdinal,Node> >(epCoordinates);
891 TEUCHOS_TEST_FOR_EXCEPT(doubleEpCoords->NumVectors() != Teuchos::as<int>(coordinates->getNumVectors()));
892 }
893 #endif
894
895 // check for Xpetra coordinates vector
896 if(paramList.isType<decltype(coordinates)>("Coordinates")) {
897 coordinates = paramList.get<decltype(coordinates)>("Coordinates");
898 }
899
900 TEUCHOS_TEST_FOR_EXCEPT(Teuchos::is_null(coordinates));
901 return coordinates;
902 }
903
904 }; // class Utilities (specialization SC=double LO=GO=int)
905
906#endif // HAVE_MUELU_EPETRA
907
908
909
939 long ExtractNonSerializableData(const Teuchos::ParameterList& inList, Teuchos::ParameterList& serialList, Teuchos::ParameterList& nonSerialList);
940
941
945 void TokenizeStringAndStripWhiteSpace(const std::string & stream, std::vector<std::string> & tokenList, const char* token = ",");
946
949 bool IsParamMuemexVariable(const std::string& name);
950
953 bool IsParamValidVariable(const std::string& name);
954
955#ifdef HAVE_MUELU_EPETRA
960 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
961 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
962 EpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Epetra_CrsMatrix>& A) {
963 typedef Xpetra::EpetraCrsMatrixT<GlobalOrdinal, Node> XECrsMatrix;
964 typedef Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrix;
965 typedef Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrixWrap;
966
967 RCP<XCrsMatrix> Atmp = rcp(new XECrsMatrix(A));
968 return rcp(new XCrsMatrixWrap(Atmp));
969 }
970
975 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
976 RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
977 EpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Epetra_MultiVector>& V) {
978 return rcp(new Xpetra::EpetraMultiVectorT<GlobalOrdinal, Node>(V));
979 }
980#endif
981
986 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
987 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
988 TpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >& Atpetra) {
989 typedef Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> XTCrsMatrix;
990 typedef Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrix;
991 typedef Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrixWrap;
992
993 RCP<XCrsMatrix> Atmp = rcp(new XTCrsMatrix(Atpetra));
994 return rcp(new XCrsMatrixWrap(Atmp));
995 }
996
1024 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1025 void leftRghtDofScalingWithinNode(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> & Amat, size_t blkSize, size_t nSweeps, Teuchos::ArrayRCP<Scalar> & rowScaling, Teuchos::ArrayRCP<Scalar> & colScaling) {
1026
1027 LocalOrdinal nBlks = (Amat.getRowMap()->getLocalNumElements())/blkSize;
1028
1029 Teuchos::ArrayRCP<Scalar> rowScaleUpdate(blkSize);
1030 Teuchos::ArrayRCP<Scalar> colScaleUpdate(blkSize);
1031
1032
1033 for (size_t i = 0; i < blkSize; i++) rowScaling[i] = 1.0;
1034 for (size_t i = 0; i < blkSize; i++) colScaling[i] = 1.0;
1035
1036 for (size_t k = 0; k < nSweeps; k++) {
1037 LocalOrdinal row = 0;
1038 for (size_t i = 0; i < blkSize; i++) rowScaleUpdate[i] = 0.0;
1039
1040 for (LocalOrdinal i = 0; i < nBlks; i++) {
1041 for (size_t j = 0; j < blkSize; j++) {
1042 Teuchos::ArrayView<const LocalOrdinal> cols;
1043 Teuchos::ArrayView<const Scalar> vals;
1044 Amat.getLocalRowView(row, cols, vals);
1045
1046 for (size_t kk = 0; kk < Teuchos::as<size_t>(vals.size()); kk++) {
1047 size_t modGuy = (cols[kk]+1)%blkSize;
1048 if (modGuy == 0) modGuy = blkSize;
1049 modGuy--;
1050 rowScaleUpdate[j] += rowScaling[j]*(Teuchos::ScalarTraits<Scalar>::magnitude(vals[kk]))*colScaling[modGuy];
1051 }
1052 row++;
1053 }
1054 }
1055 // combine information across processors
1056 Teuchos::ArrayRCP<Scalar> tempUpdate(blkSize);
1057 Teuchos::reduceAll(*(Amat.getRowMap()->getComm()), Teuchos::REDUCE_SUM, (LocalOrdinal) blkSize, rowScaleUpdate.getRawPtr(), tempUpdate.getRawPtr());
1058 for (size_t i = 0; i < blkSize; i++) rowScaleUpdate[i] = tempUpdate[i];
1059
1060 /* We want to scale by sqrt(1/rowScaleUpdate), but we'll */
1061 /* normalize things by the minimum rowScaleUpdate. That is, the */
1062 /* largest scaling is always one (as normalization is arbitrary).*/
1063
1064 Scalar minUpdate = Teuchos::ScalarTraits<Scalar>::magnitude((rowScaleUpdate[0]/rowScaling[0])/rowScaling[0]);
1065
1066 for (size_t i = 1; i < blkSize; i++) {
1067 Scalar temp = (rowScaleUpdate[i]/rowScaling[i])/rowScaling[i];
1068 if ( Teuchos::ScalarTraits<Scalar>::magnitude(temp) < Teuchos::ScalarTraits<Scalar>::magnitude(minUpdate))
1069 minUpdate = Teuchos::ScalarTraits<Scalar>::magnitude(temp);
1070 }
1071 for (size_t i = 0; i < blkSize; i++) rowScaling[i] *= sqrt(minUpdate / rowScaleUpdate[i]);
1072
1073 row = 0;
1074 for (size_t i = 0; i < blkSize; i++) colScaleUpdate[i] = 0.0;
1075
1076 for (LocalOrdinal i = 0; i < nBlks; i++) {
1077 for (size_t j = 0; j < blkSize; j++) {
1078 Teuchos::ArrayView<const LocalOrdinal> cols;
1079 Teuchos::ArrayView<const Scalar> vals;
1080 Amat.getLocalRowView(row, cols, vals);
1081 for (size_t kk = 0; kk < Teuchos::as<size_t>(vals.size()); kk++) {
1082 size_t modGuy = (cols[kk]+1)%blkSize;
1083 if (modGuy == 0) modGuy = blkSize;
1084 modGuy--;
1085 colScaleUpdate[modGuy] += colScaling[modGuy]* (Teuchos::ScalarTraits<Scalar>::magnitude(vals[kk])) *rowScaling[j];
1086 }
1087 row++;
1088 }
1089 }
1090 Teuchos::reduceAll(*(Amat.getRowMap()->getComm()), Teuchos::REDUCE_SUM, (LocalOrdinal) blkSize, colScaleUpdate.getRawPtr(), tempUpdate.getRawPtr());
1091 for (size_t i = 0; i < blkSize; i++) colScaleUpdate[i] = tempUpdate[i];
1092
1093 /* We want to scale by sqrt(1/colScaleUpdate), but we'll */
1094 /* normalize things by the minimum colScaleUpdate. That is, the */
1095 /* largest scaling is always one (as normalization is arbitrary).*/
1096
1097
1098 minUpdate = Teuchos::ScalarTraits<Scalar>::magnitude((colScaleUpdate[0]/colScaling[0])/colScaling[0]);
1099
1100 for (size_t i = 1; i < blkSize; i++) {
1101 Scalar temp = (colScaleUpdate[i]/colScaling[i])/colScaling[i];
1102 if ( Teuchos::ScalarTraits<Scalar>::magnitude(temp) < Teuchos::ScalarTraits<Scalar>::magnitude(minUpdate))
1103 minUpdate = Teuchos::ScalarTraits<Scalar>::magnitude(temp);
1104 }
1105 for (size_t i = 0; i < blkSize; i++) colScaling[i] *= sqrt(minUpdate/colScaleUpdate[i]);
1106 }
1107 }
1108
1113 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1114 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1115 TpetraFECrs_To_XpetraMatrix(const Teuchos::RCP<Tpetra::FECrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >& Atpetra) {
1116 typedef typename Tpetra::FECrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::crs_matrix_type tpetra_crs_matrix_type;
1117 typedef Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> XTCrsMatrix;
1118 typedef Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrix;
1119 typedef Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrixWrap;
1120
1121 RCP<XCrsMatrix> Atmp = rcp(new XTCrsMatrix(rcp_dynamic_cast<tpetra_crs_matrix_type>(Atpetra)));
1122 return rcp(new XCrsMatrixWrap(Atmp));
1123 }
1124
1129 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1130 RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1131 TpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >& Vtpetra) {
1132 return rcp(new Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>(Vtpetra));
1133 }
1134
1139 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1140 RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1141 TpetraFEMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Tpetra::FEMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >& Vtpetra) {
1142 typedef Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> MV;
1143 RCP<const MV> Vmv = Teuchos::rcp_dynamic_cast<const MV>(Vtpetra);
1144 return rcp(new Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>(Vmv));
1145 }
1146
1148 template<class T>
1149 std::string toString(const T& what) {
1150 std::ostringstream buf;
1151 buf << what;
1152 return buf.str();
1153 }
1154
1155#ifdef HAVE_MUELU_EPETRA
1160 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1161 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1162 EpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Epetra_CrsMatrix>& A);
1163
1168 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1169 RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1170 EpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Epetra_MultiVector>& V);
1171#endif
1172
1177 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1178 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1179 TpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >& Atpetra);
1180
1185 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1186 RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1187 TpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >& Vtpetra);
1188
1189 // Generates a communicator whose only members are other ranks of the baseComm on my node
1190 Teuchos::RCP<const Teuchos::Comm<int> > GenerateNodeComm(RCP<const Teuchos::Comm<int> > & baseComm, int &NodeId, const int reductionFactor);
1191
1192 // Lower case string
1193 std::string lowerCase (const std::string& s);
1194
1195} //namespace MueLu
1196
1197#define MUELU_UTILITIES_SHORT
1198#endif // MUELU_UTILITIES_DECL_HPP
Tpetra::KokkosCompat::KokkosSerialWrapperNode EpetraNode
MueLu::DefaultLocalOrdinal LocalOrdinal
MueLu::DefaultScalar Scalar
MueLu::DefaultGlobalOrdinal GlobalOrdinal
MueLu::DefaultNode Node
int NumMyElements() const
const Epetra_Map & RowMap() const
int ExtractMyRowView(int MyRow, int &NumEntries, double *&Values, int *&Indices) const
Exception indicating invalid cast attempted.
Exception throws to report errors in the internal logical of the program.
static const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & MV2TpetraMV(const MultiVector &vec)
static RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraBlockCrs(RCP< const Matrix > Op)
static RCP< Epetra_MultiVector > MV2NonConstEpetraMV(RCP< MultiVector > vec)
static RCP< Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraRow(RCP< Operator > Op)
Teuchos::ScalarTraits< Scalar >::magnitudeType Magnitude
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV2(MultiVector &vec)
static void MyOldScaleMatrix(Matrix &Op, const Teuchos::ArrayRCP< const Scalar > &scalingVector, bool doInverse=true, bool doFillComplete=true, bool doOptimizeStorage=true)
static Epetra_MultiVector & MV2NonConstEpetraMV(MultiVector &vec)
static void MyOldScaleMatrix_Epetra(Matrix &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool, bool)
static RCP< Matrix > Transpose(Matrix &Op, bool=false, const std::string &label=std::string(), const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Transpose a Xpetra::Matrix.
static void MyOldScaleMatrix_Tpetra(Matrix &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV(RCP< MultiVector > vec)
static RCP< const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraRow(RCP< const Operator > Op)
static Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2NonConstTpetraBlockCrs(Matrix &Op)
static RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraCrs(RCP< Matrix > Op)
static const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2TpetraBlockCrs(const Matrix &Op)
static RCP< const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2TpetraMV(RCP< MultiVector > const vec)
Helper utility to pull out the underlying Tpetra objects from an Xpetra object.
static const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2TpetraCrs(const Matrix &Op)
static RCP< const Epetra_MultiVector > MV2EpetraMV(RCP< MultiVector > const vec)
Helper utility to pull out the underlying Epetra objects from an Xpetra object.
static const Epetra_Map & Map2EpetraMap(const Map &map)
static const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > Map2TpetraMap(const Map &map)
static RCP< const Epetra_CrsMatrix > Op2EpetraCrs(RCP< const Matrix > Op)
static const Epetra_MultiVector & MV2EpetraMV(const MultiVector &vec)
static Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2NonConstTpetraCrs(Matrix &Op)
static RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraCrs(RCP< const Matrix > Op)
Xpetra::EpetraMultiVectorT< GlobalOrdinal, Node > EpetraMultiVector
static RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::magnitudeType, LocalOrdinal, GlobalOrdinal, Node > > ExtractCoordinatesFromParameterList(ParameterList &paramList)
Extract coordinates from parameter list and return them in a Xpetra::MultiVector.
static RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::magnitudeType, LocalOrdinal, GlobalOrdinal, Node > > RealValuedToScalarMultiVector(RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::coordinateType, LocalOrdinal, GlobalOrdinal, Node > > X)
static RCP< Epetra_CrsMatrix > Op2NonConstEpetraCrs(RCP< Matrix > Op)
static Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & MV2NonConstTpetraMV(MultiVector &vec)
static const Epetra_CrsMatrix & Op2EpetraCrs(const Matrix &Op)
static RCP< Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraBlockCrs(RCP< Matrix > Op)
MueLu utility class.
static const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2TpetraCrs(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op)
static RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraCrs(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraBlockCrs(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraRow(RCP< Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
Teuchos::ScalarTraits< Scalar >::magnitudeType Magnitude
static const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & MV2TpetraMV(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
static void MyOldScaleMatrix_Tpetra(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
static Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2NonConstTpetraCrs(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op)
static RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Transpose(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, bool optimizeTranspose=false, const std::string &label=std::string(), const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
static const Epetra_MultiVector & MV2EpetraMV(const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
static RCP< const Epetra_CrsMatrix > Op2EpetraCrs(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV2(Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
static Epetra_CrsMatrix & Op2NonConstEpetraCrs(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op)
static const Epetra_Map & Map2EpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
static Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2NonConstTpetraBlockCrs(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op)
static RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > RealValuedToScalarMultiVector(RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::coordinateType, LocalOrdinal, GlobalOrdinal, Node > > X)
static const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > Map2TpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
static Epetra_MultiVector & MV2NonConstEpetraMV(Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
static RCP< const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2TpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > const vec)
Helper utility to pull out the underlying Tpetra objects from an Xpetra object.
static void MyOldScaleMatrix(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< const Scalar > &scalingVector, bool doInverse=true, bool doFillComplete=true, bool doOptimizeStorage=true)
static Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & MV2NonConstTpetraMV(Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
static void MyOldScaleMatrix_Epetra(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec)
static RCP< Epetra_MultiVector > MV2NonConstEpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec)
static RCP< Epetra_CrsMatrix > Op2NonConstEpetraCrs(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraCrs(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraBlockCrs(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraRow(RCP< const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static const Epetra_CrsMatrix & Op2EpetraCrs(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op)
static RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::magnitudeType, LocalOrdinal, GlobalOrdinal, Node > > ExtractCoordinatesFromParameterList(ParameterList &paramList)
static RCP< const Epetra_MultiVector > MV2EpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > const vec)
Helper utility to pull out the underlying Epetra objects from an Xpetra object.
static const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2TpetraBlockCrs(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op)
Namespace for MueLu classes and methods.
bool IsParamMuemexVariable(const std::string &name)
long ExtractNonSerializableData(const Teuchos::ParameterList &inList, Teuchos::ParameterList &serialList, Teuchos::ParameterList &nonSerialList)
Extract non-serializable data from level-specific sublists and move it to a separate parameter list.
std::string lowerCase(const std::string &s)
RCP< Xpetra::MultiVector< SC, LO, GO, NO > > TpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP< Tpetra::MultiVector< SC, LO, GO, NO > > &Vtpetra)
void leftRghtDofScalingWithinNode(const Xpetra::Matrix< SC, LO, GO, NO > &Atpetra, size_t blkSize, size_t nSweeps, Teuchos::ArrayRCP< SC > &rowScaling, Teuchos::ArrayRCP< SC > &colScaling)
RCP< Xpetra::MultiVector< SC, LO, GO, NO > > TpetraFEMultiVector_To_XpetraMultiVector(const Teuchos::RCP< Tpetra::FEMultiVector< SC, LO, GO, NO > > &Vtpetra)
RCP< Xpetra::MultiVector< SC, LO, GO, NO > > EpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP< Epetra_MultiVector > &V)
Tpetra::KokkosClassic::DefaultNode::DefaultNodeType DefaultNode
RCP< Xpetra::Matrix< SC, LO, GO, NO > > EpetraCrs_To_XpetraMatrix(const Teuchos::RCP< Epetra_CrsMatrix > &A)
bool IsParamValidVariable(const std::string &name)
void TokenizeStringAndStripWhiteSpace(const std::string &stream, std::vector< std::string > &tokenList, const char *delimChars)
RCP< Xpetra::Matrix< SC, LO, GO, NO > > TpetraCrs_To_XpetraMatrix(const Teuchos::RCP< Tpetra::CrsMatrix< SC, LO, GO, NO > > &Atpetra)
RCP< Xpetra::Matrix< SC, LO, GO, NO > > TpetraFECrs_To_XpetraMatrix(const Teuchos::RCP< Tpetra::FECrsMatrix< SC, LO, GO, NO > > &Atpetra)
RCP< Xpetra::CrsMatrixWrap< SC, LO, GO, NO > > Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap(RCP< Epetra_CrsMatrix > &epAB)
std::string toString(const T &what)
Little helper function to convert non-string types to strings.
Teuchos::RCP< const Teuchos::Comm< int > > GenerateNodeComm(RCP< const Teuchos::Comm< int > > &baseComm, int &NodeId, const int reductionFactor)