MueLu  Version of the Day
MueLu_Utilities_def.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_DEF_HPP
47 #define MUELU_UTILITIES_DEF_HPP
48 
49 #include <Teuchos_DefaultComm.hpp>
50 #include <Teuchos_ParameterList.hpp>
51 
52 #include "MueLu_ConfigDefs.hpp"
53 
54 #ifdef HAVE_MUELU_EPETRA
55 # ifdef HAVE_MPI
56 # include "Epetra_MpiComm.h"
57 # endif
58 #endif
59 
60 #if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT)
61 #include <EpetraExt_MatrixMatrix.h>
62 #include <EpetraExt_RowMatrixOut.h>
64 #include <EpetraExt_CrsMatrixIn.h>
66 #include <EpetraExt_BlockMapIn.h>
67 #include <Xpetra_EpetraUtils.hpp>
68 #include <Xpetra_EpetraMultiVector.hpp>
69 #include <EpetraExt_BlockMapOut.h>
70 #endif
71 
72 #ifdef HAVE_MUELU_TPETRA
73 #include <MatrixMarket_Tpetra.hpp>
74 #include <Tpetra_RowMatrixTransposer.hpp>
75 #include <TpetraExt_MatrixMatrix.hpp>
76 #include <Xpetra_TpetraMultiVector.hpp>
77 #include <Xpetra_TpetraCrsMatrix.hpp>
78 #include <Xpetra_TpetraBlockCrsMatrix.hpp>
79 #endif
80 
81 #ifdef HAVE_MUELU_EPETRA
82 #include <Xpetra_EpetraMap.hpp>
83 #endif
84 
85 #include <Xpetra_BlockedCrsMatrix.hpp>
86 //#include <Xpetra_DefaultPlatform.hpp>
87 #include <Xpetra_IO.hpp>
88 #include <Xpetra_Import.hpp>
89 #include <Xpetra_ImportFactory.hpp>
90 #include <Xpetra_Map.hpp>
91 #include <Xpetra_MapFactory.hpp>
92 #include <Xpetra_Matrix.hpp>
93 #include <Xpetra_MatrixFactory.hpp>
94 #include <Xpetra_MultiVector.hpp>
95 #include <Xpetra_MultiVectorFactory.hpp>
96 #include <Xpetra_Operator.hpp>
97 #include <Xpetra_Vector.hpp>
98 #include <Xpetra_VectorFactory.hpp>
99 
100 #include <Xpetra_MatrixMatrix.hpp>
101 
102 #include <MueLu_Utilities_decl.hpp>
103 #if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_ML)
104 #include <ml_operator.h>
105 #include <ml_epetra_utils.h>
106 #endif
107 
108 namespace MueLu {
109 
110 #ifdef HAVE_MUELU_EPETRA
111  //using Xpetra::EpetraCrsMatrix; // TODO: mv in Xpetra_UseShortNamesScalar
112  //using Xpetra::EpetraMultiVector;
113 #endif
114 
115 #ifdef HAVE_MUELU_EPETRA
116  template<typename SC,typename LO,typename GO,typename NO>
117  RCP<Xpetra::CrsMatrixWrap<SC,LO,GO,NO> > Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap(RCP<Epetra_CrsMatrix> &epAB){
118  return Xpetra::Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap<SC,LO,GO,NO>(epAB);
119  }
120 #endif
121 
122 #ifdef HAVE_MUELU_EPETRA
123  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
124  RCP<const Epetra_MultiVector> Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MV2EpetraMV(const RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > vec) {
125  RCP<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> >(vec);
126  if (tmpVec == Teuchos::null)
127  throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
128  return tmpVec->getEpetra_MultiVector();
129  }
130 
131  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
132  RCP<Epetra_MultiVector> Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MV2NonConstEpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > vec) {
133  RCP<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> >(vec);
134  if (tmpVec == Teuchos::null)
135  throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
136  return tmpVec->getEpetra_MultiVector();
137  }
138 
139  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
140  Epetra_MultiVector& Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MV2NonConstEpetraMV(Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> &vec) {
141  const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> & tmpVec = dynamic_cast<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> &>(vec);
142  return *(tmpVec.getEpetra_MultiVector());
143  }
144 
145  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
146  const Epetra_MultiVector& Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MV2EpetraMV(const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec) {
147  const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> & tmpVec = dynamic_cast<const Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> &>(vec);
148  return *(tmpVec.getEpetra_MultiVector());
149  }
150 
151  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
152  RCP<const Epetra_CrsMatrix> Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Op2EpetraCrs(RCP<const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op) {
153  RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
154  if (crsOp == Teuchos::null)
155  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
156  const RCP<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>>& tmp_ECrsMtx = rcp_dynamic_cast<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
157  if (tmp_ECrsMtx == Teuchos::null)
158  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
159  return tmp_ECrsMtx->getEpetra_CrsMatrix();
160  }
161 
162  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
163  RCP<Epetra_CrsMatrix> Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Op2NonConstEpetraCrs(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op) {
164  RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
165  if (crsOp == Teuchos::null)
166  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
167  const RCP<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
168  if (tmp_ECrsMtx == Teuchos::null)
169  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
170  return tmp_ECrsMtx->getEpetra_CrsMatrixNonConst();
171  }
172 
173  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
174  const Epetra_CrsMatrix& Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Op2EpetraCrs(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op) {
175  try {
176  const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>& crsOp = dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(Op);
177  try {
178  const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>& tmp_ECrsMtx = dynamic_cast<const Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>&>(*crsOp.getCrsMatrix());
179  return *tmp_ECrsMtx.getEpetra_CrsMatrix();
180  } catch (std::bad_cast&) {
181  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
182  }
183  } catch (std::bad_cast&) {
184  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
185  }
186  }
187 
188  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
189  Epetra_CrsMatrix& Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Op2NonConstEpetraCrs(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op) {
190  try {
191  Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>& crsOp = dynamic_cast<Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(Op);
192  try {
193  Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>& tmp_ECrsMtx = dynamic_cast<Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>&>(*crsOp.getCrsMatrix());
194  return *tmp_ECrsMtx.getEpetra_CrsMatrixNonConst();
195  } catch (std::bad_cast&) {
196  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
197  }
198  } catch (std::bad_cast&) {
199  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
200  }
201  }
202 
203  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
204  const Epetra_Map& Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Map2EpetraMap(const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node>& map) {
205  RCP<const Xpetra::EpetraMapT<GlobalOrdinal,Node> > xeMap = rcp_dynamic_cast<const Xpetra::EpetraMapT<GlobalOrdinal,Node> >(rcpFromRef(map));
206  if (xeMap == Teuchos::null)
207  throw Exceptions::BadCast("Utilities::Map2EpetraMap : Cast from Xpetra::Map to Xpetra::EpetraMap failed");
208  return xeMap->getEpetra_Map();
209  }
210 #endif
211 
212 #ifdef HAVE_MUELU_TPETRA
213  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
214  RCP<const Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
215  Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MV2TpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > const vec) {
216  RCP<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(vec);
217  if (tmpVec == Teuchos::null)
218  throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
219  return tmpVec->getTpetra_MultiVector();
220  }
221 
222  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
223  RCP<Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> > Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MV2NonConstTpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > vec) {
224  RCP<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmpVec = rcp_dynamic_cast<Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(vec);
225  if (tmpVec == Teuchos::null)
226  throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
227  return tmpVec->getTpetra_MultiVector();
228  }
229 
230  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
231  Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> & Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MV2NonConstTpetraMV(Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec) {
232  const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec = dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(vec);
233  return *(tmpVec.getTpetra_MultiVector());
234  }
235 
236  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
237  RCP<Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> > Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MV2NonConstTpetraMV2(Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> &vec) {
238  const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec = dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(vec);
239  return tmpVec.getTpetra_MultiVector();
240  }
241 
242  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
243  const Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>&
244  Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MV2TpetraMV(const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec) {
245  const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec = dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(vec);
246  return *(tmpVec.getTpetra_MultiVector());
247  }
248 
249  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
250  RCP<const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Op2TpetraCrs(RCP<const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > Op) {
251  // Get the underlying Tpetra Mtx
252  RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
253  if (crsOp == Teuchos::null)
254  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
255  const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
256  if (tmp_ECrsMtx == Teuchos::null)
257  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
258  return tmp_ECrsMtx->getTpetra_CrsMatrix();
259  }
260 
261  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
262  RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Op2NonConstTpetraCrs(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op) {
263  RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
264  if (crsOp == Teuchos::null)
265  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
266  const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &tmp_ECrsMtx = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsOp->getCrsMatrix());
267  if (tmp_ECrsMtx == Teuchos::null)
268  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
269  return tmp_ECrsMtx->getTpetra_CrsMatrixNonConst();
270  }
271 
272  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
273  const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Op2TpetraCrs(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op) {
274  try {
275  const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>& crsOp = dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(Op);
276  try {
277  const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx = dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(*crsOp.getCrsMatrix());
278  return *tmp_ECrsMtx.getTpetra_CrsMatrix();
279  } catch (std::bad_cast&) {
280  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
281  }
282  } catch (std::bad_cast&) {
283  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
284  }
285  }
286 
287  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
288  Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Op2NonConstTpetraCrs(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op) {
289  try {
290  Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>& crsOp = dynamic_cast<Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(Op);
291  try {
292  Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx = dynamic_cast<Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(*crsOp.getCrsMatrix());
293  return *tmp_ECrsMtx.getTpetra_CrsMatrixNonConst();
294  } catch (std::bad_cast&) {
295  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
296  }
297  } catch (std::bad_cast&) {
298  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
299  }
300  }
301 
302  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
303  RCP<const Tpetra::RowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Op2TpetraRow(RCP<const Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > Op) {
304  RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
305  if (crsOp == Teuchos::null)
306  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
307 
308  RCP<const Xpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsMat = crsOp->getCrsMatrix();
309  const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_Crs = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
310  RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_BlockCrs;
311  if(!tmp_Crs.is_null()) {
312  return tmp_Crs->getTpetra_CrsMatrixNonConst();
313  }
314  else {
315  tmp_BlockCrs= rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
316  if (tmp_BlockCrs.is_null())
317  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
318  return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
319  }
320  }
321 
322  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
323  RCP<Tpetra::RowMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Op2NonConstTpetraRow(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op) {
324  RCP<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsOp = rcp_dynamic_cast<const Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(Op);
325  if (crsOp == Teuchos::null)
326  throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
327 
328  RCP<const Xpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > crsMat = crsOp->getCrsMatrix();
329  const RCP<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_Crs = rcp_dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
330  RCP<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > tmp_BlockCrs;
331  if(!tmp_Crs.is_null()) {
332  return tmp_Crs->getTpetra_CrsMatrixNonConst();
333  }
334  else {
335  tmp_BlockCrs= rcp_dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(crsMat);
336  if (tmp_BlockCrs.is_null())
337  throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
338  return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
339  }
340  }
341 
342 
343  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
344  const RCP<const Tpetra::Map<LocalOrdinal, GlobalOrdinal,Node> > Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Map2TpetraMap(const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node>& map) {
345  const RCP<const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node>>& tmp_TMap = rcp_dynamic_cast<const Xpetra::TpetraMap<LocalOrdinal,GlobalOrdinal,Node> >(rcpFromRef(map));
346  if (tmp_TMap == Teuchos::null)
347  throw Exceptions::BadCast("Utilities::Map2TpetraMap : Cast from Xpetra::Map to Xpetra::TpetraMap failed");
348  return tmp_TMap->getTpetra_Map();
349  }
350 #endif
351 
352  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
353  void Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MyOldScaleMatrix(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Teuchos::ArrayRCP<const Scalar>& scalingVector, bool doInverse,
354  bool doFillComplete,
355  bool doOptimizeStorage)
356  {
357  Scalar one = Teuchos::ScalarTraits<Scalar>::one();
358  Teuchos::ArrayRCP<Scalar> sv(scalingVector.size());
359  if (doInverse) {
360  for (int i = 0; i < scalingVector.size(); ++i)
361  sv[i] = one / scalingVector[i];
362  } else {
363  for (int i = 0; i < scalingVector.size(); ++i)
364  sv[i] = scalingVector[i];
365  }
366 
367  switch (Op.getRowMap()->lib()) {
368  case Xpetra::UseTpetra:
369  MyOldScaleMatrix_Tpetra(Op, sv, doFillComplete, doOptimizeStorage);
370  break;
371 
372  case Xpetra::UseEpetra:
373  MyOldScaleMatrix_Epetra(Op, sv, doFillComplete, doOptimizeStorage);
374  break;
375 
376  default:
377  throw Exceptions::RuntimeError("Only Epetra and Tpetra matrices can be scaled.");
378  }
379  }
380 
381  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
382  void Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MyOldScaleMatrix_Epetra(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& /* Op */, const Teuchos::ArrayRCP<Scalar>& /* scalingVector */, bool /* doFillComplete */, bool /* doOptimizeStorage */) {
383  throw Exceptions::RuntimeError("MyOldScaleMatrix_Epetra: Epetra needs SC=double and LO=GO=int.");
384  }
385 
386  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
387  void Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MyOldScaleMatrix_Tpetra(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Teuchos::ArrayRCP<Scalar>& scalingVector,
388  bool doFillComplete,
389  bool doOptimizeStorage)
390  {
391 #ifdef HAVE_MUELU_TPETRA
392  try {
393  Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& tpOp = Op2NonConstTpetraCrs(Op);
394 
395  const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rowMap = tpOp.getRowMap();
396  const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > domainMap = tpOp.getDomainMap();
397  const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > rangeMap = tpOp.getRangeMap();
398 
399  size_t maxRowSize = tpOp.getNodeMaxNumRowEntries();
400  if (maxRowSize == Teuchos::as<size_t>(-1)) // hasn't been determined yet
401  maxRowSize = 20;
402 
403 
404  if (tpOp.isFillComplete())
405  tpOp.resumeFill();
406 
407  if (Op.isLocallyIndexed() == true) {
408  typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::local_inds_host_view_type cols;
409  typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::values_host_view_type vals;
410 
411  for (size_t i = 0; i < rowMap->getNodeNumElements(); ++i) {
412  tpOp.getLocalRowView(i, cols, vals);
413 
414  size_t nnz = tpOp.getNumEntriesInLocalRow(i);
415  typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::nonconst_values_host_view_type scaledVals("ScaledVals",nnz);
416 
417  for (size_t j = 0; j < nnz; ++j) {
418  scaledVals[j] = scalingVector[i]*vals[j];
419  }
420 
421  if (nnz > 0) {
422  tpOp.replaceLocalValues(i, cols, scaledVals);
423  }
424  } //for (size_t i=0; ...
425 
426  } else {
427  typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::global_inds_host_view_type cols;
428  typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::values_host_view_type vals;
429 
430  for (size_t i = 0; i < rowMap->getNodeNumElements(); ++i) {
431  GlobalOrdinal gid = rowMap->getGlobalElement(i);
432  tpOp.getGlobalRowView(gid, cols, vals);
433  size_t nnz = tpOp.getNumEntriesInGlobalRow(gid);
434  typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::nonconst_values_host_view_type scaledVals("ScaledVals",nnz);
435 
436  // FIXME FIXME FIXME FIXME FIXME FIXME
437  for (size_t j = 0; j < nnz; ++j)
438  scaledVals[j] = scalingVector[i]*vals[j]; //FIXME i or gid?
439 
440  if (nnz > 0) {
441  tpOp.replaceGlobalValues(gid, cols, scaledVals);
442  }
443  } //for (size_t i=0; ...
444  }
445 
446  if (doFillComplete) {
447  if (domainMap == Teuchos::null || rangeMap == Teuchos::null)
448  throw Exceptions::RuntimeError("In Utilities::Scaling: cannot fillComplete because the domain and/or range map hasn't been defined");
449 
450  RCP<Teuchos::ParameterList> params = rcp(new Teuchos::ParameterList());
451  params->set("Optimize Storage", doOptimizeStorage);
452  params->set("No Nonlocal Changes", true);
453  Op.fillComplete(Op.getDomainMap(), Op.getRangeMap(), params);
454  }
455  } catch(...) {
456  throw Exceptions::RuntimeError("Only Tpetra::CrsMatrix types can be scaled (Err.1)");
457  }
458 #else
459  throw Exceptions::RuntimeError("Matrix scaling is not possible because Tpetra has not been enabled.");
460 #endif
461  } //MyOldScaleMatrix_Tpetra()
462 
463  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
464  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
466  Transpose (Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, bool /* optimizeTranspose */,const std::string & label,const Teuchos::RCP<Teuchos::ParameterList> &params) {
467 #if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT)
468  std::string TorE = "epetra";
469 #else
470  std::string TorE = "tpetra";
471 #endif
472 
473 #if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT)
474  try {
476  (void) epetraOp; // silence "unused variable" compiler warning
477  } catch (...) {
478  TorE = "tpetra";
479  }
480 #endif
481 
482 #ifdef HAVE_MUELU_TPETRA
483  if (TorE == "tpetra") {
484  try {
485  const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& tpetraOp = Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Op2TpetraCrs(Op);
486 
487  RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > A;
488  Tpetra::RowMatrixTransposer<Scalar, LocalOrdinal, GlobalOrdinal, Node> transposer(rcpFromRef(tpetraOp),label); //more than meets the eye
489 
490  {
491  using Teuchos::ParameterList;
492  using Teuchos::rcp;
493  RCP<ParameterList> transposeParams = params.is_null () ?
494  rcp (new ParameterList) :
495  rcp (new ParameterList (*params));
496  transposeParams->set ("sort", false);
497  A = transposer.createTranspose (transposeParams);
498  }
499 
500  RCP<Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AA = rcp(new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(A) );
501  RCP<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AAA = rcp_implicit_cast<Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >(AA);
502  RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AAAA = rcp( new Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node>(AAA) );
503  if (!AAAA->isFillComplete())
504  AAAA->fillComplete(Op.getRangeMap(), Op.getDomainMap());
505 
506  if (Op.IsView("stridedMaps"))
507  AAAA->CreateView("stridedMaps", Teuchos::rcpFromRef(Op), true/*doTranspose*/);
508 
509  return AAAA;
510 
511  } catch (std::exception& e) {
512  std::cout << "threw exception '" << e.what() << "'" << std::endl;
513  throw Exceptions::RuntimeError("Utilities::Transpose failed, perhaps because matrix is not a Crs matrix");
514  }
515  } //if
516 #endif
517 
518  // Epetra case
519  std::cout << "Utilities::Transpose() not implemented for Epetra" << std::endl;
520  return Teuchos::null;
521 
522  } // Transpose
523 
524 
525  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
526  RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
528  RealValuedToScalarMultiVector(RCP<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LocalOrdinal,GlobalOrdinal,Node> > X) {
529  RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Xscalar;
530 #if defined(HAVE_XPETRA_TPETRA) && (defined(HAVE_TPETRA_INST_COMPLEX_DOUBLE) || defined(HAVE_TPETRA_INST_COMPLEX_FLOAT))
531  typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType real_type;
532  // Need to cast the real-valued multivector to Scalar=complex
533  if ((typeid(Scalar).name() == typeid(std::complex<double>).name()) ||
534  (typeid(Scalar).name() == typeid(std::complex<float>).name())) {
535  Xscalar = Xpetra::MultiVectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(X->getMap(),X->getNumVectors());
536  size_t numVecs = X->getNumVectors();
537  for (size_t j=0;j<numVecs;j++) {
538  Teuchos::ArrayRCP<const real_type> XVec = X->getData(j);
539  Teuchos::ArrayRCP<Scalar> XVecScalar = Xscalar->getDataNonConst(j);
540  for(size_t i = 0; i < static_cast<size_t>(XVec.size()); ++i)
541  XVecScalar[i]=XVec[i];
542  }
543  } else
544 #endif
545  Xscalar = rcp_dynamic_cast<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(X);
546  return Xscalar;
547  }
548 
549 
550  template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
551  RCP<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LocalOrdinal,GlobalOrdinal,Node> >
553  ExtractCoordinatesFromParameterList (ParameterList& paramList) {
554 
555  RCP<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LocalOrdinal,GlobalOrdinal,Node> > coordinates = Teuchos::null;
556 
557  // check whether coordinates are contained in parameter list
558  if(paramList.isParameter ("Coordinates") == false)
559  return coordinates;
560 
561 #if defined(HAVE_MUELU_TPETRA)
562  // only Tpetra code
563 
564  // define Tpetra::MultiVector type with Scalar=float only if
565  // * ETI is turned off, since then the compiler will instantiate it automatically OR
566  // * Tpetra is instantiated on Scalar=float
567 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT)
568  typedef Tpetra::MultiVector<float, LocalOrdinal, GlobalOrdinal, Node> tfMV;
569  RCP<tfMV> floatCoords = Teuchos::null;
570 #endif
571 
572  // define Tpetra::MultiVector type with Scalar=double only if
573  // * ETI is turned off, since then the compiler will instantiate it automatically OR
574  // * Tpetra is instantiated on Scalar=double
575 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_DOUBLE)
576  typedef Tpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType, LocalOrdinal, GlobalOrdinal, Node> tdMV;
577  RCP<tdMV> doubleCoords = Teuchos::null;
578  if (paramList.isType<RCP<tdMV> >("Coordinates")) {
579  // Coordinates are stored as a double vector
580  doubleCoords = paramList.get<RCP<tdMV> >("Coordinates");
581  paramList.remove("Coordinates");
582  }
583 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT)
584  else if (paramList.isType<RCP<tfMV> >("Coordinates")) {
585  // check if coordinates are stored as a float vector
586  floatCoords = paramList.get<RCP<tfMV> >("Coordinates");
587  paramList.remove("Coordinates");
588  doubleCoords = rcp(new tdMV(floatCoords->getMap(), floatCoords->getNumVectors()));
589  deep_copy(*doubleCoords, *floatCoords);
590  }
591 #endif
592  // We have the coordinates in a Tpetra double vector
593  if(doubleCoords != Teuchos::null) {
594  //rcp(new Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>(Vtpetra));
595  coordinates = Teuchos::rcp_dynamic_cast<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LocalOrdinal,GlobalOrdinal,Node> >(MueLu::TpetraMultiVector_To_XpetraMultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LocalOrdinal,GlobalOrdinal,Node>(doubleCoords));
596  TEUCHOS_TEST_FOR_EXCEPT(Teuchos::is_null(coordinates));
597  TEUCHOS_TEST_FOR_EXCEPT(doubleCoords->getNumVectors() != coordinates->getNumVectors());
598  }
599 #else
600  // coordinates usually are stored as double vector
601  // Tpetra is not instantiated on scalar=double
602  throw Exceptions::RuntimeError("ExtractCoordinatesFromParameterList: The coordinates vector in parameter list is expected to be a Tpetra multivector with SC=double or float.");
603 #endif
604 #endif // endif HAVE_TPETRA
605 
606  // check for Xpetra coordinates vector
607  if(paramList.isType<decltype(coordinates)>("Coordinates")) {
608  coordinates = paramList.get<decltype(coordinates)>("Coordinates");
609  }
610 
611  return coordinates;
612  } // ExtractCoordinatesFromParameterList
613 
614 } //namespace MueLu
615 
616 #define MUELU_UTILITIES_SHORT
617 #endif // MUELU_UTILITIES_DEF_HPP
618 
619 // LocalWords: LocalOrdinal
Exception indicating invalid cast attempted.
static const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > Map2TpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
MueLu::DefaultLocalOrdinal LocalOrdinal
static RCP< const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2TpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > const vec)
Helper utility to pull out the underlying Tpetra objects from an Xpetra object.
static void MyOldScaleMatrix(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< const Scalar > &scalingVector, bool doInverse=true, bool doFillComplete=true, bool doOptimizeStorage=true)
static RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraCrs(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static const Epetra_Map & Map2EpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
static RCP< Epetra_CrsMatrix > Op2NonConstEpetraCrs(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::magnitudeType, LocalOrdinal, GlobalOrdinal, Node > > ExtractCoordinatesFromParameterList(ParameterList &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.
Namespace for MueLu classes and methods.
MueLu::DefaultNode Node
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec)
MueLu::DefaultScalar Scalar
MueLu::DefaultGlobalOrdinal GlobalOrdinal
static RCP< Epetra_MultiVector > MV2NonConstEpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec)
static RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > RealValuedToScalarMultiVector(RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::magnitudeType, LocalOrdinal, GlobalOrdinal, Node > > X)
static RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraCrs(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
void deep_copy(const DynRankView< DT, DP... > &dst, typename ViewTraits< DT, DP... >::const_value_type &value, typename std::enable_if< std::is_same< typename ViewTraits< DT, DP... >::specialize, void >::value >::type *=nullptr)
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV2(Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &vec)
static RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Transpose(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, bool optimizeTranspose=false, const std::string &label=std::string(), const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
static RCP< const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraRow(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
RCP< Xpetra::CrsMatrixWrap< SC, LO, GO, NO > > Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap(RCP< Epetra_CrsMatrix > &epAB)
static void MyOldScaleMatrix_Epetra(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
static RCP< const Epetra_CrsMatrix > Op2EpetraCrs(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
Exception throws to report errors in the internal logical of the program.
static RCP< Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraRow(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static void MyOldScaleMatrix_Tpetra(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)