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#ifdef HAVE_MUELU_TPETRA
58#include <Xpetra_TpetraBlockCrsMatrix.hpp>
59#include <Xpetra_TpetraOperator.hpp>
60#endif
61#include <Xpetra_BlockedCrsMatrix_fwd.hpp>
62#include <Xpetra_CrsGraph_fwd.hpp>
63#include <Xpetra_CrsGraphFactory_fwd.hpp>
64#include <Xpetra_CrsMatrix_fwd.hpp>
65#include <Xpetra_CrsMatrixWrap_fwd.hpp>
66#include <Xpetra_Map_fwd.hpp>
67#include <Xpetra_MapFactory_fwd.hpp>
68#include <Xpetra_Matrix_fwd.hpp>
69#include <Xpetra_MatrixFactory_fwd.hpp>
70#include <Xpetra_MultiVector_fwd.hpp>
71#include <Xpetra_MultiVectorFactory_fwd.hpp>
72#include <Xpetra_Operator_fwd.hpp>
73#include <Xpetra_Vector_fwd.hpp>
74#include <Xpetra_VectorFactory_fwd.hpp>
75#include <Xpetra_ExportFactory.hpp>
76
77#include <Xpetra_Import.hpp>
78#include <Xpetra_ImportFactory.hpp>
79#include <Xpetra_MatrixMatrix.hpp>
80
81#ifdef HAVE_MUELU_EPETRA
82#include <Xpetra_EpetraCrsMatrix_fwd.hpp>
83
84// needed because of inlined function
85//TODO: remove inline function?
86#include <Xpetra_EpetraCrsMatrix.hpp>
87#include <Xpetra_CrsMatrixWrap.hpp>
88
89#endif
90
91#include "MueLu_Exceptions.hpp"
92
93#ifdef HAVE_MUELU_EPETRAEXT
96class Epetra_Vector;
98#endif
99
100#ifdef HAVE_MUELU_TPETRA
101#include <Tpetra_CrsMatrix.hpp>
102#include <Tpetra_BlockCrsMatrix.hpp>
103#include <Tpetra_BlockCrsMatrix_Helpers.hpp>
104#include <Tpetra_FECrsMatrix.hpp>
105#include <Tpetra_RowMatrixTransposer.hpp>
106#include <Tpetra_Map.hpp>
107#include <Tpetra_MultiVector.hpp>
108#include <Tpetra_FEMultiVector.hpp>
109#include <Xpetra_TpetraCrsMatrix_fwd.hpp>
110#include <Xpetra_TpetraMultiVector_fwd.hpp>
111#endif
112
113#include <MueLu_UtilitiesBase.hpp>
114
115
116namespace MueLu {
117
118#ifdef HAVE_MUELU_EPETRA
119 //defined after Utilities class
120 template<typename SC,typename LO,typename GO,typename NO>
121 RCP<Xpetra::CrsMatrixWrap<SC,LO,GO,NO> >
122 Convert_Epetra_CrsMatrix_ToXpetra_CrsMatrixWrap(RCP<Epetra_CrsMatrix> &epAB);
123
124 template<typename SC,typename LO,typename GO,typename NO>
125 RCP<Xpetra::Matrix<SC, LO, GO, NO> >
126 EpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Epetra_CrsMatrix>& A);
127
128 template<typename SC,typename LO,typename GO,typename NO>
129 RCP<Xpetra::MultiVector<SC, LO, GO, NO> >
130 EpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Epetra_MultiVector>& V);
131#endif
132
133#ifdef HAVE_MUELU_TPETRA
134 template<typename SC,typename LO,typename GO,typename NO>
135 RCP<Xpetra::Matrix<SC, LO, GO, NO> >
136 TpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Tpetra::CrsMatrix<SC, LO, GO, NO> >& Atpetra);
137
138 template<typename SC,typename LO,typename GO,typename NO>
139 RCP<Xpetra::Matrix<SC, LO, GO, NO> >
140 TpetraFECrs_To_XpetraMatrix(const Teuchos::RCP<Tpetra::FECrsMatrix<SC, LO, GO, NO> >& Atpetra);
141
142 template<typename SC,typename LO,typename GO,typename NO>
143 RCP<Xpetra::MultiVector<SC, LO, GO, NO> >
144 TpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Tpetra::MultiVector<SC, LO, GO, NO> >& Vtpetra);
145
146 template<typename SC,typename LO,typename GO,typename NO>
147 RCP<Xpetra::MultiVector<SC, LO, GO, NO> >
148 TpetraFEMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Tpetra::FEMultiVector<SC, LO, GO, NO> >& Vtpetra);
149
150 template<typename SC,typename LO,typename GO,typename NO>
151 void leftRghtDofScalingWithinNode(const Xpetra::Matrix<SC,LO,GO,NO> & Atpetra, size_t blkSize, size_t nSweeps, Teuchos::ArrayRCP<SC> & rowScaling, Teuchos::ArrayRCP<SC> & colScaling);
152#endif
153
161 template <class Scalar,
164 class Node = DefaultNode>
165 class Utilities : public UtilitiesBase<Scalar, LocalOrdinal, GlobalOrdinal, Node> {
166#undef MUELU_UTILITIES_SHORT
167 //#include "MueLu_UseShortNames.hpp"
168
169 public:
170 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType Magnitude;
171
172#ifdef HAVE_MUELU_EPETRA
174 // @{
175 static RCP<const Epetra_MultiVector> MV2EpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > const vec);
176 static RCP< Epetra_MultiVector> MV2NonConstEpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > vec);
177
178 static const Epetra_MultiVector& MV2EpetraMV(const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
179 static Epetra_MultiVector& MV2NonConstEpetraMV(Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
180
181 static RCP<const Epetra_CrsMatrix> Op2EpetraCrs(RCP<const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
182 static RCP< Epetra_CrsMatrix> Op2NonConstEpetraCrs(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
183
184 static const Epetra_CrsMatrix& Op2EpetraCrs(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
185 static Epetra_CrsMatrix& Op2NonConstEpetraCrs(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
186
187 static const Epetra_Map& Map2EpetraMap(const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node>& map);
188 // @}
189#endif
190
191#ifdef HAVE_MUELU_TPETRA
193 static RCP<const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2TpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > const vec);
194 static RCP< Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2NonConstTpetraMV(RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > vec);
195 static RCP< Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > MV2NonConstTpetraMV2(Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
196
197 static const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& MV2TpetraMV(const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
198 static Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& MV2NonConstTpetraMV(Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& vec);
199
200 static RCP<const Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2TpetraCrs(RCP<const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
201 static RCP< Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2NonConstTpetraCrs(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
202
203 static const Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2TpetraCrs(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
204 static Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2NonConstTpetraCrs(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
205
206 static RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2TpetraBlockCrs(RCP<const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
207 static RCP< Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2NonConstTpetraBlockCrs(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
208
209 static const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2TpetraBlockCrs(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
210 static Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2NonConstTpetraBlockCrs(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op);
211
212
213
214 static RCP<const Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2TpetraRow(RCP<const Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
215 static RCP< Tpetra::RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op2NonConstTpetraRow(RCP<Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op);
216
217
218 static const RCP<const Tpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > Map2TpetraMap(const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node>& map);
219#endif
220
221 static RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Crs2Op(RCP<Xpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Op) { return UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Crs2Op(Op); }
224 static Teuchos::ArrayRCP<Scalar> GetMatrixDiagonal(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetMatrixDiagonal(A); }
225 static RCP<Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > GetMatrixDiagonalInverse(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, Magnitude tol = Teuchos::ScalarTraits<Scalar>::eps()*100, Scalar valReplacement = Teuchos::ScalarTraits<Scalar>::zero()) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetMatrixDiagonalInverse(A,tol,valReplacement); }
226 static Teuchos::RCP<Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > GetLumpedMatrixDiagonal(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> const &A, const bool doReciprocal=false, Magnitude tol = Teuchos::ScalarTraits<Scalar>::eps()*100, Scalar tolReplacement = Teuchos::ScalarTraits<Scalar>::zero(), const bool replaceSingleEntryRowWithZero = false, const bool useAverageAbsDiagVal = false)
229 static Teuchos::ArrayRCP<Magnitude> GetMatrixMaxMinusOffDiagonal(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetMatrixMaxMinusOffDiagonal(A); }
230 static Teuchos::ArrayRCP<Magnitude> GetMatrixMaxMinusOffDiagonal(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, const Xpetra::Vector<LocalOrdinal,LocalOrdinal,GlobalOrdinal,Node> &BlockNumber) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetMatrixMaxMinusOffDiagonal(A,BlockNumber); }
231
232 static Teuchos::RCP<Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > GetInverse(Teuchos::RCP<const Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > v, Magnitude tol = Teuchos::ScalarTraits<Scalar>::eps()*100, Scalar tolReplacement = Teuchos::ScalarTraits<Scalar>::zero()) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetInverse(v,tol,tolReplacement); }
233 static Teuchos::Array<Magnitude> ResidualNorm(const Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& X, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& RHS) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::ResidualNorm(Op,X,RHS); }
234 static Teuchos::Array<Magnitude> ResidualNorm(const Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& X, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& RHS, Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Resid) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::ResidualNorm(Op,X,RHS,Resid); }
235 static RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Residual(const Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& X, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& RHS) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Residual(Op,X,RHS); }
236 static void Residual(const Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& X, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& RHS, Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Resid) { MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Residual(Op,X,RHS,Resid);}
238 static typename Teuchos::ScalarTraits<Scalar>::magnitudeType Distance2(const Teuchos::Array<Teuchos::ArrayRCP<const Scalar>>& v, LocalOrdinal i0, LocalOrdinal i1) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Distance2(v,i0,i1); }
239 static typename Teuchos::ScalarTraits<Scalar>::magnitudeType Distance2(const Teuchos::ArrayView<double> & weight,const Teuchos::Array<Teuchos::ArrayRCP<const Scalar>>& v, LocalOrdinal i0, LocalOrdinal i1) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Distance2(weight,v,i0,i1); }
240 static Teuchos::ArrayRCP<const bool> DetectDirichletRows(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, const Magnitude& tol = Teuchos::ScalarTraits<Scalar>::magnitude(0.), const bool count_twos_as_dirichlet=false) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::DetectDirichletRows(A,tol,count_twos_as_dirichlet); }
241 static Teuchos::ArrayRCP<const bool> DetectDirichletRowsExt(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, bool & bHasZeroDiagonal, const Magnitude& tol = Teuchos::ScalarTraits<Scalar>::zero()) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::DetectDirichletRowsExt(A,bHasZeroDiagonal,tol); }
242 static void FindNonZeros(const Teuchos::ArrayRCP<const Scalar> vals,Teuchos::ArrayRCP<bool> nonzeros) {return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::FindNonZeros(vals,nonzeros);}
243 static void DetectDirichletColsAndDomains(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, const Teuchos::ArrayRCP<bool>& dirichletRows, Teuchos::ArrayRCP<bool> dirichletCols, Teuchos::ArrayRCP<bool> dirichletDomains) {MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::DetectDirichletColsAndDomains(A,dirichletRows,dirichletCols,dirichletDomains); }
244
245 static void ApplyRowSumCriterion(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, const Magnitude rowSumTol, Teuchos::ArrayRCP<bool>& dirichletRows) {return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::ApplyRowSumCriterion(A,rowSumTol,dirichletRows); }
246 static void ApplyRowSumCriterion(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, const Xpetra::Vector<LocalOrdinal,LocalOrdinal,GlobalOrdinal,Node> &BlockNumber,const Magnitude rowSumTol, Teuchos::ArrayRCP<bool>& dirichletRows) {return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::ApplyRowSumCriterion(A,BlockNumber,rowSumTol,dirichletRows); }
247
249
250 static Scalar PowerMethod(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, bool scaleByDiag = true,
251 LocalOrdinal niters = 10, Magnitude tolerance = 1e-2, bool verbose = false, unsigned int seed = 123) {
253 }
254
255 static Scalar PowerMethod(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, const RCP<Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node>> &invDiag,
256 LocalOrdinal niters = 10, Magnitude tolerance = 1e-2, bool verbose = false, unsigned int seed = 123) {
258 }
259
260 static Scalar Frobenius(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& B) {
262 }
263
264 static void MyOldScaleMatrix(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Teuchos::ArrayRCP<const Scalar>& scalingVector, bool doInverse = true,
265 bool doFillComplete = true, bool doOptimizeStorage = true);
266
267 static void MyOldScaleMatrix_Epetra(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Teuchos::ArrayRCP<Scalar>& scalingVector,
269 static void MyOldScaleMatrix_Tpetra(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Teuchos::ArrayRCP<Scalar>& scalingVector,
271
272 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);
273
274 static RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > RealValuedToScalarMultiVector(RCP<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::coordinateType,LocalOrdinal,GlobalOrdinal,Node> > X);
275
276 static RCP<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LocalOrdinal,GlobalOrdinal,Node> > ExtractCoordinatesFromParameterList(ParameterList& paramList);
277
278 static void FindDirichletRows(Teuchos::RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > & A, std::vector<LocalOrdinal>& dirichletRows,bool count_twos_as_dirichlet=false) {
280 }
281
282
283 static void ApplyOAZToMatrixRows(Teuchos::RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& A,const std::vector<LocalOrdinal>& dirichletRows) {
285 }
286
287 static void ApplyOAZToMatrixRows(Teuchos::RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& A,const Teuchos::ArrayRCP<const bool>& dirichletRows) {
289 }
290
291 static void ZeroDirichletRows(Teuchos::RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& A,const std::vector<LocalOrdinal>& dirichletRows, Scalar replaceWith=Teuchos::ScalarTraits<Scalar>::zero()) {
293 }
294
295 static void ZeroDirichletRows(Teuchos::RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& A,const Teuchos::ArrayRCP<const bool>& dirichletRows, Scalar replaceWith=Teuchos::ScalarTraits<Scalar>::zero()) {
297 }
298
299 static void ZeroDirichletRows(Teuchos::RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& X,const Teuchos::ArrayRCP<const bool>& dirichletRows,Scalar replaceWith=Teuchos::ScalarTraits<Scalar>::zero()) {
301 }
302
303 static void ZeroDirichletCols(Teuchos::RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& A,const Teuchos::ArrayRCP<const bool>& dirichletCols, Scalar replaceWith=Teuchos::ScalarTraits<Scalar>::zero()) {
305 }
306
307 }; // class Utilities
308
310
311#ifdef HAVE_MUELU_EPETRA
321 template <>
322 class Utilities<double,int,int,Xpetra::EpetraNode> : public UtilitiesBase<double,int,int,Xpetra::EpetraNode> {
323 public:
324 typedef double Scalar;
325 typedef int LocalOrdinal;
326 typedef int GlobalOrdinal;
327 typedef Xpetra::EpetraNode Node;
328 typedef Teuchos::ScalarTraits<Scalar>::magnitudeType Magnitude;
329
330 private:
331 using CrsGraph = Xpetra::CrsGraph<LocalOrdinal,GlobalOrdinal,Node>;
332 using CrsMatrixWrap = Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
333 using CrsMatrix = Xpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
334 using Matrix = Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
335 using Operator = Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
336 using Vector = Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
337 using MultiVector = Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
338 using Map = Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node>;
339#ifdef HAVE_MUELU_TPETRA
340 using TpetraOperator = Xpetra::TpetraOperator<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
341#endif
342#ifdef HAVE_MUELU_EPETRA
343 using EpetraMap = Xpetra::EpetraMapT<GlobalOrdinal,Node>;
344 using EpetraMultiVector = Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node>;
345 using EpetraCrsMatrix = Xpetra::EpetraCrsMatrixT<GlobalOrdinal,Node>;
346#endif
347 public:
348
349#ifdef HAVE_MUELU_EPETRA
351 // @{
354 if (tmpVec == Teuchos::null)
355 throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
356 return tmpVec->getEpetra_MultiVector();
357 }
360 if (tmpVec == Teuchos::null)
361 throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::EpetraMultiVector failed");
362 return tmpVec->getEpetra_MultiVector();
363 }
364
366 const EpetraMultiVector& tmpVec = dynamic_cast<const EpetraMultiVector&>(vec);
367 return *(tmpVec.getEpetra_MultiVector());
368 }
370 const EpetraMultiVector& tmpVec = dynamic_cast<const EpetraMultiVector&>(vec);
371 return *(tmpVec.getEpetra_MultiVector());
372 }
373
376 if (crsOp == Teuchos::null)
377 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
379 if (tmp_ECrsMtx == Teuchos::null)
380 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
381 return tmp_ECrsMtx->getEpetra_CrsMatrix();
382 }
385 if (crsOp == Teuchos::null)
386 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
388 if (tmp_ECrsMtx == Teuchos::null)
389 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
390 return tmp_ECrsMtx->getEpetra_CrsMatrixNonConst();
391 }
392
393 static const Epetra_CrsMatrix& Op2EpetraCrs(const Matrix& Op) {
394 try {
395 const CrsMatrixWrap& crsOp = dynamic_cast<const CrsMatrixWrap&>(Op);
396 try {
397 const EpetraCrsMatrix& tmp_ECrsMtx = dynamic_cast<const EpetraCrsMatrix&>(*crsOp.getCrsMatrix());
398 return *tmp_ECrsMtx.getEpetra_CrsMatrix();
399 } catch (std::bad_cast&) {
400 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
401 }
402 } catch (std::bad_cast&) {
403 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
404 }
405 }
407 try {
408 CrsMatrixWrap& crsOp = dynamic_cast<CrsMatrixWrap&>(Op);
409 try {
410 EpetraCrsMatrix& tmp_ECrsMtx = dynamic_cast<EpetraCrsMatrix&>(*crsOp.getCrsMatrix());
411 return *tmp_ECrsMtx.getEpetra_CrsMatrixNonConst();
412 } catch (std::bad_cast&) {
413 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::EpetraCrsMatrix failed");
414 }
415 } catch (std::bad_cast&) {
416 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
417 }
418 }
419
420 static const Epetra_Map& Map2EpetraMap(const Map& map) {
422 if (xeMap == Teuchos::null)
423 throw Exceptions::BadCast("Utilities::Map2EpetraMap : Cast from Xpetra::Map to Xpetra::EpetraMap failed");
424 return xeMap->getEpetra_Map();
425 }
426 // @}
427#endif
428
429#ifdef HAVE_MUELU_TPETRA
432#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
433 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
434 throw Exceptions::RuntimeError("MV2TpetraMV: Tpetra has not been compiled with support for LO=GO=int.");
435#else
437 if (tmpVec == Teuchos::null)
438 throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
439 return tmpVec->getTpetra_MultiVector();
440#endif
441 }
443#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
444 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
445 throw Exceptions::RuntimeError("MV2NonConstTpetraMV: Tpetra has not been compiled with support for LO=GO=int.");
446#else
448 if (tmpVec == Teuchos::null)
449 throw Exceptions::BadCast("Cast from Xpetra::MultiVector to Xpetra::TpetraMultiVector failed");
450 return tmpVec->getTpetra_MultiVector();
451#endif
452
453 }
455#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
456 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
457 throw Exceptions::RuntimeError("MV2NonConstTpetraMV2: Tpetra has not been compiled with support for LO=GO=int.");
458#else
459 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec = dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(vec);
460 return tmpVec.getTpetra_MultiVector();
461#endif
462 }
463
464 static const Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& MV2TpetraMV(const MultiVector& vec) {
465#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
466 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
467 throw Exceptions::RuntimeError("MV2TpetraMV: Tpetra has not been compiled with support for LO=GO=int.");
468#else
469 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec = dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(vec);
470 return *(tmpVec.getTpetra_MultiVector());
471#endif
472 }
473 static Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& MV2NonConstTpetraMV(MultiVector& vec) {
474#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
475 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
476 throw Exceptions::RuntimeError("MV2NonConstTpetraMV: Tpetra has not been compiled with support for LO=GO=int.");
477#else
478 const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmpVec = dynamic_cast<const Xpetra::TpetraMultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(vec);
479 return *(tmpVec.getTpetra_MultiVector());
480#endif
481 }
482
484#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
485 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
486 throw Exceptions::RuntimeError("Op2TpetraCrs: Tpetra has not been compiled with support for LO=GO=int.");
487#else
488 // Get the underlying Tpetra Mtx
490 if (crsOp == Teuchos::null)
491 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
493 if (tmp_ECrsMtx == Teuchos::null)
494 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
495 return tmp_ECrsMtx->getTpetra_CrsMatrix();
496#endif
497 }
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
504 if (crsOp == Teuchos::null)
505 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
507 if (tmp_ECrsMtx == Teuchos::null)
508 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
509 return tmp_ECrsMtx->getTpetra_CrsMatrixNonConst();
510#endif
511 };
512
513 static const Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2TpetraCrs(const Matrix& Op) {
514#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
515 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
516 throw Exceptions::RuntimeError("Op2TpetraCrs: Tpetra has not been compiled with support for LO=GO=int.");
517#else
518 try {
519 const CrsMatrixWrap& crsOp = dynamic_cast<const CrsMatrixWrap&>(Op);
520 try {
521 const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx = dynamic_cast<const Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(*crsOp.getCrsMatrix());
522 return *tmp_ECrsMtx.getTpetra_CrsMatrix();
523 } catch (std::bad_cast&) {
524 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
525 }
526 } catch (std::bad_cast&) {
527 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
528 }
529#endif
530 }
531 static Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2NonConstTpetraCrs(Matrix& Op) {
532#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
533 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
534 throw Exceptions::RuntimeError("Op2NonConstTpetraCrs: Tpetra has not been compiled with support for LO=GO=int.");
535#else
536 try {
537 CrsMatrixWrap& crsOp = dynamic_cast<CrsMatrixWrap&>(Op);
538 try {
539 Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx = dynamic_cast<Xpetra::TpetraCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(*crsOp.getCrsMatrix());
540 return *tmp_ECrsMtx.getTpetra_CrsMatrixNonConst();
541 } catch (std::bad_cast&) {
542 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix failed");
543 }
544 } catch (std::bad_cast&) {
545 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
546 }
547#endif
548 }
549
550
552#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
553 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
554 throw Exceptions::RuntimeError("Op2TpetraBlockCrs: Tpetra has not been compiled with support for LO=GO=int.");
555#else
556 // Get the underlying Tpetra Mtx
558 if (crsOp == Teuchos::null)
559 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
561 if (tmp_ECrsMtx == Teuchos::null)
562 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
563 return tmp_ECrsMtx->getTpetra_BlockCrsMatrix();
564#endif
565 }
566
568#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
569 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
570 throw Exceptions::RuntimeError("Op2NonConstTpetraBlockCrs: Tpetra has not been compiled with support for LO=GO=int.");
571#else
573 if (crsOp == Teuchos::null)
574 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
576 if (tmp_ECrsMtx == Teuchos::null)
577 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
578 return tmp_ECrsMtx->getTpetra_BlockCrsMatrixNonConst();
579#endif
580 };
581
582 static const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2TpetraBlockCrs(const Matrix& Op) {
583#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
584 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
585 throw Exceptions::RuntimeError("Op2TpetraBlockCrs: Tpetra has not been compiled with support for LO=GO=int.");
586#else
587 try {
588 const CrsMatrixWrap& crsOp = dynamic_cast<const CrsMatrixWrap&>(Op);
589 try {
590 const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx = dynamic_cast<const Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(*crsOp.getCrsMatrix());
591 return *tmp_ECrsMtx.getTpetra_BlockCrsMatrix();
592 } catch (std::bad_cast&) {
593 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
594 }
595 } catch (std::bad_cast&) {
596 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
597 }
598#endif
599 }
600 static Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op2NonConstTpetraBlockCrs(Matrix& Op) {
601#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
602 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
603 throw Exceptions::RuntimeError("Op2NonConstTpetraCrs: Tpetra has not been compiled with support for LO=GO=int.");
604#else
605 try {
606 CrsMatrixWrap& crsOp = dynamic_cast<CrsMatrixWrap&>(Op);
607 try {
608 Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& tmp_ECrsMtx = dynamic_cast<Xpetra::TpetraBlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>&>(*crsOp.getCrsMatrix());
609 return *tmp_ECrsMtx.getTpetra_BlockCrsMatrixNonConst();
610 } catch (std::bad_cast&) {
611 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraBlockCrsMatrix failed");
612 }
613 } catch (std::bad_cast&) {
614 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
615 }
616#endif
617 }
618
619
621#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
622 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
623 throw Exceptions::RuntimeError("Op2TpetraRow: Tpetra has not been compiled with support for LO=GO=int.");
624#else
626 if (!mat.is_null()) {
628 if (crsOp == Teuchos::null)
629 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
630
631 RCP<const CrsMatrix> crsMat = crsOp->getCrsMatrix();
634 if(!tmp_Crs.is_null()) {
635 return tmp_Crs->getTpetra_CrsMatrixNonConst();
636 }
637 else {
639 if (tmp_BlockCrs.is_null())
640 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
641 return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
642 }
643 } else {
647 return tRow;
648 }
649#endif
650 }
651
653#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
654 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
655 throw Exceptions::RuntimeError("Op2NonConstTpetraRow: Tpetra has not been compiled with support for LO=GO=int.");
656#else
658 if (!mat.is_null()) {
660 if (crsOp == Teuchos::null)
661 throw Exceptions::BadCast("Cast from Xpetra::Matrix to Xpetra::CrsMatrixWrap failed");
662
663 RCP<const CrsMatrix> crsMat = crsOp->getCrsMatrix();
666 if(!tmp_Crs.is_null()) {
667 return tmp_Crs->getTpetra_CrsMatrixNonConst();
668 }
669 else {
671 if (tmp_BlockCrs.is_null())
672 throw Exceptions::BadCast("Cast from Xpetra::CrsMatrix to Xpetra::TpetraCrsMatrix and Xpetra::TpetraBlockCrsMatrix failed");
673 return tmp_BlockCrs->getTpetra_BlockCrsMatrixNonConst();
674 }
675 } else {
679 return tRow;
680 }
681#endif
682 };
683
684
686#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
687 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
688 throw Exceptions::RuntimeError("Map2TpetraMap: Tpetra has not been compiled with support for LO=GO=int.");
689#else
691 if (tmp_TMap == Teuchos::null)
692 throw Exceptions::BadCast("Utilities::Map2TpetraMap : Cast from Xpetra::Map to Xpetra::TpetraMap failed");
693 return tmp_TMap->getTpetra_Map();
694#endif
695 };
696#endif
697
702 static RCP<Vector> GetMatrixDiagonalInverse(const Matrix& A, Magnitude tol = Teuchos::ScalarTraits<Scalar>::eps()*100, Scalar valReplacement = Teuchos::ScalarTraits<Scalar>::zero()) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetMatrixDiagonalInverse(A,tol,valReplacement); }
703 static Teuchos::RCP<Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > GetLumpedMatrixDiagonal(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> const &A, const bool doReciprocal=false, Magnitude tol = Teuchos::ScalarTraits<Scalar>::eps()*100, Scalar tolReplacement = Teuchos::ScalarTraits<Scalar>::zero(), const bool replaceSingleEntryRowWithZero = false, const bool useAverageAbsDiagVal = false)
706 static Teuchos::ArrayRCP<Magnitude> GetMatrixMaxMinusOffDiagonal(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetMatrixMaxMinusOffDiagonal(A); }
707 static Teuchos::ArrayRCP<Magnitude> GetMatrixMaxMinusOffDiagonal(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, const Xpetra::Vector<LocalOrdinal,LocalOrdinal,GlobalOrdinal,Node> &BlockNumber) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetMatrixMaxMinusOffDiagonal(A,BlockNumber); }
708 static RCP<Vector> GetInverse(Teuchos::RCP<const Vector> v, Magnitude tol = Teuchos::ScalarTraits<Scalar>::eps()*100, Scalar tolReplacement = Teuchos::ScalarTraits<Scalar>::zero()) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::GetInverse(v,tol,tolReplacement); }
709 static Teuchos::Array<Magnitude> ResidualNorm(const Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& X, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& RHS) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::ResidualNorm(Op,X,RHS); }
710 static Teuchos::Array<Magnitude> ResidualNorm(const Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& X, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& RHS, Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Resid) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::ResidualNorm(Op,X,RHS,Resid); }
711 static RCP<MultiVector> Residual(const Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& X, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& RHS) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Residual(Op,X,RHS); }
712 static void Residual(const Xpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& X, const Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& RHS, Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Resid) { MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Residual(Op,X,RHS,Resid);}
714 static Teuchos::ScalarTraits<Scalar>::magnitudeType Distance2(const Teuchos::Array<Teuchos::ArrayRCP<const Scalar>>& v, LocalOrdinal i0, LocalOrdinal i1) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Distance2(v,i0,i1); }
715 static Teuchos::ScalarTraits<Scalar>::magnitudeType Distance2(const Teuchos::ArrayView<double> &weight, const Teuchos::Array<Teuchos::ArrayRCP<const Scalar>>& v, LocalOrdinal i0, LocalOrdinal i1) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Distance2(weight,v,i0,i1); }
716 static Teuchos::ArrayRCP<const bool> DetectDirichletRows(const Matrix& A, const Magnitude& tol = Teuchos::ScalarTraits<Scalar>::zero(), const bool count_twos_as_dirichlet=false) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::DetectDirichletRows(A,tol,count_twos_as_dirichlet); }
717 static Teuchos::ArrayRCP<const bool> DetectDirichletRowsExt(const Matrix& A, bool & bHasZeroDiagonal, const Magnitude& tol = Teuchos::ScalarTraits<Scalar>::zero()) { return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::DetectDirichletRowsExt(A,bHasZeroDiagonal,tol); }
718 static void FindNonZeros(const Teuchos::ArrayRCP<const Scalar> vals,Teuchos::ArrayRCP<bool> nonzeros) {return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::FindNonZeros(vals,nonzeros);}
719 static void DetectDirichletColsAndDomains(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, const Teuchos::ArrayRCP<bool>& dirichletRows, Teuchos::ArrayRCP<bool> dirichletCols, Teuchos::ArrayRCP<bool> dirichletDomains) {MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::DetectDirichletColsAndDomains(A,dirichletRows,dirichletCols,dirichletDomains); }
720
721
722 static void ApplyRowSumCriterion(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, const Magnitude rowSumTol, Teuchos::ArrayRCP<bool>& dirichletRows) {return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::ApplyRowSumCriterion(A,rowSumTol,dirichletRows); }
723 static void ApplyRowSumCriterion(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, const Xpetra::Vector<LocalOrdinal,LocalOrdinal,GlobalOrdinal,Node> &BlockNumber,const Magnitude rowSumTol, Teuchos::ArrayRCP<bool>& dirichletRows) {return MueLu::UtilitiesBase<Scalar,LocalOrdinal,GlobalOrdinal,Node>::ApplyRowSumCriterion(A,BlockNumber,rowSumTol,dirichletRows); }
725
726 static Scalar PowerMethod(const Matrix& A, bool scaleByDiag = true,
727 LocalOrdinal niters = 10, Magnitude tolerance = 1e-2, bool verbose = false, unsigned int seed = 123) {
729 }
730
735
736 static Scalar Frobenius(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A, const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& B) {
738 }
739
740 static void MyOldScaleMatrix(Matrix& Op, const Teuchos::ArrayRCP<const Scalar>& scalingVector, bool doInverse = true,
741 bool doFillComplete = true, bool doOptimizeStorage = true) {
742 Scalar one = Teuchos::ScalarTraits<Scalar>::one();
743 Teuchos::ArrayRCP<Scalar> sv(scalingVector.size());
744 if (doInverse) {
745 for (int i = 0; i < scalingVector.size(); ++i)
746 sv[i] = one / scalingVector[i];
747 } else {
748 for (int i = 0; i < scalingVector.size(); ++i)
749 sv[i] = scalingVector[i];
750 }
751
752 switch (Op.getRowMap()->lib()) {
753 case Xpetra::UseTpetra:
755 break;
756
757 case Xpetra::UseEpetra:
759 break;
760
761 default:
762 throw Exceptions::RuntimeError("Only Epetra and Tpetra matrices can be scaled.");
763 }
764 }
765
766 // TODO This is the <double,int,int> specialization
767 static void MyOldScaleMatrix_Tpetra(Matrix& Op, const Teuchos::ArrayRCP<Scalar>& scalingVector,
769#ifdef HAVE_MUELU_TPETRA
770#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
771 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
772 throw Exceptions::RuntimeError("Matrix scaling is not possible because Tpetra has not been compiled with support for LO=GO=int.");
773#else
774 try {
775 Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& tpOp = Op2NonConstTpetraCrs(Op);
776
778 const RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > domainMap = tpOp.getDomainMap();
780
781 size_t maxRowSize = tpOp.getLocalMaxNumRowEntries();
782 if (maxRowSize == Teuchos::as<size_t>(-1)) // hasn't been determined yet
783 maxRowSize = 20;
784
785 std::vector<Scalar> scaledVals(maxRowSize);
786 if (tpOp.isFillComplete())
787 tpOp.resumeFill();
788
789 if (Op.isLocallyIndexed() == true) {
790 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::local_inds_host_view_type cols;
791 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::values_host_view_type vals;
792 for (size_t i = 0; i < rowMap->getLocalNumElements(); ++i) {
793 tpOp.getLocalRowView(i, cols, vals);
794 size_t nnz = tpOp.getNumEntriesInLocalRow(i);
795 if (nnz > maxRowSize) {
796 maxRowSize = nnz;
797 scaledVals.resize(maxRowSize);
798 }
799 for (size_t j = 0; j < nnz; ++j)
800 scaledVals[j] = vals[j]*scalingVector[i];
801
802 if (nnz > 0) {
803 Teuchos::ArrayView<const LocalOrdinal> cols_view(cols.data(), nnz);
804 Teuchos::ArrayView<const Scalar> valview(&scaledVals[0], nnz);
805 tpOp.replaceLocalValues(i, cols_view, valview);
806 }
807 } //for (size_t i=0; ...
808
809 } else {
810 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::local_inds_host_view_type cols;
811 typename Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::values_host_view_type vals;
812
813 for (size_t i = 0; i < rowMap->getLocalNumElements(); ++i) {
814 GlobalOrdinal gid = rowMap->getGlobalElement(i);
815 tpOp.getGlobalRowView(gid, cols, vals);
816 size_t nnz = tpOp.getNumEntriesInGlobalRow(gid);
817 if (nnz > maxRowSize) {
818 maxRowSize = nnz;
819 scaledVals.resize(maxRowSize);
820 }
821 // FIXME FIXME FIXME FIXME FIXME FIXME
822 for (size_t j = 0; j < nnz; ++j)
823 scaledVals[j] = vals[j]*scalingVector[i]; //FIXME i or gid?
824
825 if (nnz > 0) {
826 Teuchos::ArrayView<const LocalOrdinal> cols_view(cols.data(), nnz);
827 Teuchos::ArrayView<const Scalar> valview(&scaledVals[0], nnz);
828 tpOp.replaceGlobalValues(gid, cols_view, valview);
829 }
830 } //for (size_t i=0; ...
831 }
832
833 if (doFillComplete) {
834 if (domainMap == Teuchos::null || rangeMap == Teuchos::null)
835 throw Exceptions::RuntimeError("In Utilities::Scaling: cannot fillComplete because the domain and/or range map hasn't been defined");
836
837 RCP<Teuchos::ParameterList> params = rcp(new Teuchos::ParameterList());
838 params->set("Optimize Storage", doOptimizeStorage);
839 params->set("No Nonlocal Changes", true);
840 Op.fillComplete(Op.getDomainMap(), Op.getRangeMap(), params);
841 }
842 } catch(...) {
843 throw Exceptions::RuntimeError("Only Tpetra::CrsMatrix types can be scaled (Err.1)");
844 }
845#endif
846#else
847 throw Exceptions::RuntimeError("Matrix scaling is not possible because Tpetra has not been enabled.");
848#endif
849 }
850
851 static void MyOldScaleMatrix_Epetra (Matrix& Op, const Teuchos::ArrayRCP<Scalar>& scalingVector, bool /* doFillComplete */, bool /* doOptimizeStorage */) {
852#ifdef HAVE_MUELU_EPETRA
853 try {
854 //const Epetra_CrsMatrix& epOp = Utilities<double,int,int>::Op2NonConstEpetraCrs(Op);
856
857 Epetra_Map const &rowMap = epOp.RowMap();
858 int nnz;
859 double *vals;
860 int *cols;
861
862 for (int i = 0; i < rowMap.NumMyElements(); ++i) {
863 epOp.ExtractMyRowView(i, nnz, vals, cols);
864 for (int j = 0; j < nnz; ++j)
865 vals[j] *= scalingVector[i];
866 }
867
868 } catch (...){
869 throw Exceptions::RuntimeError("Only Epetra_CrsMatrix types can be scaled");
870 }
871#else
872 throw Exceptions::RuntimeError("Matrix scaling is not possible because Epetra has not been enabled.");
873#endif // HAVE_MUELU_EPETRA
874 }
875
881 static RCP<Matrix> Transpose(Matrix& Op, bool /* optimizeTranspose */ = false,const std::string & label = std::string(),const Teuchos::RCP<Teuchos::ParameterList> &params=Teuchos::null) {
882 switch (Op.getRowMap()->lib()) {
883 case Xpetra::UseTpetra: {
884#ifdef HAVE_MUELU_TPETRA
885#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
886 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
887 throw Exceptions::RuntimeError("Utilities::Transpose: Tpetra is not compiled with LO=GO=int. Add TPETRA_INST_INT_INT:BOOL=ON to your configuration!");
888#else
889 using Helpers = Xpetra::Helpers<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
890 /***************************************************************/
891 if(Helpers::isTpetraCrs(Op)) {
892 const Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>& tpetraOp = Utilities<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Op2TpetraCrs(Op);
893
894 // Compute the transpose A of the Tpetra matrix tpetraOp.
896 Tpetra::RowMatrixTransposer<Scalar, LocalOrdinal, GlobalOrdinal, Node> transposer(rcpFromRef(tpetraOp),label);
897
898 {
899 using Teuchos::ParameterList;
900 using Teuchos::rcp;
902 rcp (new ParameterList) :
903 rcp (new ParameterList (*params));
904 transposeParams->set ("sort", false);
905 A = transposer.createTranspose(transposeParams);
906 }
907
908 RCP<Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AA = rcp(new Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(A));
910 RCP<Matrix> AAAA = rcp( new CrsMatrixWrap(AAA));
911
912 if (Op.IsView("stridedMaps"))
913 AAAA->CreateView("stridedMaps", Teuchos::rcpFromRef(Op), true/*doTranspose*/);
914
915 return AAAA;
916 }
917 /***************************************************************/
918 else if(Helpers::isTpetraBlockCrs(Op)) {
919 using BCRS = Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
920 // using CRS = Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>;
923 {
924 Tpetra::BlockCrsMatrixTransposer<Scalar, LocalOrdinal, GlobalOrdinal, Node> transposer(rcpFromRef(tpetraOp),label);
925
926 using Teuchos::ParameterList;
927 using Teuchos::rcp;
929 rcp (new ParameterList) :
930 rcp (new ParameterList (*params));
931 transposeParams->set ("sort", false);
932 At = transposer.createTranspose(transposeParams);
933 }
934
935 RCP<Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > AA = rcp(new Xpetra::TpetraBlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(At));
937 RCP<Matrix> AAAA = rcp( new CrsMatrixWrap(AAA));
938
939 if (Op.IsView("stridedMaps"))
940 AAAA->CreateView("stridedMaps", Teuchos::rcpFromRef(Op), true/*doTranspose*/);
941
942 return AAAA;
943
944 }
945 /***************************************************************/
946 else {
947 throw Exceptions::RuntimeError("Utilities::Transpose failed, perhaps because matrix is not a Crs or BlockCrs matrix");
948 }
949#endif
950#else
951 throw Exceptions::RuntimeError("Utilities::Transpose: Tpetra is not compiled!");
952#endif
953 }
954 case Xpetra::UseEpetra:
955 {
956#if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT)
957 Teuchos::TimeMonitor tm(*Teuchos::TimeMonitor::getNewTimer("ZZ Entire Transpose"));
958 // Epetra case
962 transposer.ReleaseTranspose(); // So we can keep A in Muelu...
963
967 RCP<Matrix> AAAA = rcp( new CrsMatrixWrap(AAA));
968
969 if (Op.IsView("stridedMaps"))
970 AAAA->CreateView("stridedMaps", Teuchos::rcpFromRef(Op), true/*doTranspose*/);
971
972 return AAAA;
973#else
974 throw Exceptions::RuntimeError("Epetra (Err. 2)");
975#endif
976 }
977 default:
978 throw Exceptions::RuntimeError("Only Epetra and Tpetra matrices can be transposed.");
979 }
980
981 TEUCHOS_UNREACHABLE_RETURN(Teuchos::null);
982 }
983
984 static RCP<Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LocalOrdinal,GlobalOrdinal,Node> >
989
994
995 // check whether coordinates are contained in parameter list
996 if(paramList.isParameter ("Coordinates") == false)
997 return coordinates;
998
999 #if defined(HAVE_MUELU_TPETRA)
1000 #if ( defined(EPETRA_HAVE_OMP) && defined(HAVE_TPETRA_INST_OPENMP) && defined(HAVE_TPETRA_INST_INT_INT)) || \
1001 (!defined(EPETRA_HAVE_OMP) && defined(HAVE_TPETRA_INST_SERIAL) && defined(HAVE_TPETRA_INST_INT_INT))
1002
1003 // define Tpetra::MultiVector type with Scalar=float only if
1004 // * ETI is turned off, since then the compiler will instantiate it automatically OR
1005 // * Tpetra is instantiated on Scalar=float
1006 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT)
1007 typedef Tpetra::MultiVector<float, LocalOrdinal, GlobalOrdinal, Node> tfMV;
1008 RCP<tfMV> floatCoords = Teuchos::null;
1009 #endif
1010
1011 // define Tpetra::MultiVector type with Scalar=double only if
1012 // * ETI is turned off, since then the compiler will instantiate it automatically OR
1013 // * Tpetra is instantiated on Scalar=double
1014 typedef Tpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType, LocalOrdinal, GlobalOrdinal, Node> tdMV;
1015 RCP<tdMV> doubleCoords = Teuchos::null;
1016 if (paramList.isType<RCP<tdMV> >("Coordinates")) {
1017 // Coordinates are stored as a double vector
1018 doubleCoords = paramList.get<RCP<tdMV> >("Coordinates");
1019 paramList.remove("Coordinates");
1020 }
1021 #if !defined(HAVE_TPETRA_EXPLICIT_INSTANTIATION) || defined(HAVE_TPETRA_INST_FLOAT)
1022 else if (paramList.isType<RCP<tfMV> >("Coordinates")) {
1023 // check if coordinates are stored as a float vector
1024 floatCoords = paramList.get<RCP<tfMV> >("Coordinates");
1025 paramList.remove("Coordinates");
1026 doubleCoords = rcp(new tdMV(floatCoords->getMap(), floatCoords->getNumVectors()));
1028 }
1029 #endif
1030 // We have the coordinates in a Tpetra double vector
1031 if(doubleCoords != Teuchos::null) {
1032 coordinates = Teuchos::rcp(new Xpetra::TpetraMultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType, LocalOrdinal, GlobalOrdinal, Node>(doubleCoords));
1033 TEUCHOS_TEST_FOR_EXCEPT(doubleCoords->getNumVectors() != coordinates->getNumVectors());
1034 }
1035 #endif // Tpetra instantiated on GO=int and EpetraNode
1036 #endif // endif HAVE_TPETRA
1037
1038 #if defined(HAVE_MUELU_EPETRA)
1040 if (paramList.isType<RCP<Epetra_MultiVector> >("Coordinates")) {
1041 doubleEpCoords = paramList.get<RCP<Epetra_MultiVector> >("Coordinates");
1042 paramList.remove("Coordinates");
1043 RCP<Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node> > epCoordinates = Teuchos::rcp(new Xpetra::EpetraMultiVectorT<GlobalOrdinal,Node>(doubleEpCoords));
1045 TEUCHOS_TEST_FOR_EXCEPT(doubleEpCoords->NumVectors() != Teuchos::as<int>(coordinates->getNumVectors()));
1046 }
1047 #endif
1048
1049 // check for Xpetra coordinates vector
1050 if(paramList.isType<decltype(coordinates)>("Coordinates")) {
1051 coordinates = paramList.get<decltype(coordinates)>("Coordinates");
1052 }
1053
1054 TEUCHOS_TEST_FOR_EXCEPT(Teuchos::is_null(coordinates));
1055 return coordinates;
1056 }
1057
1058 }; // class Utilities (specialization SC=double LO=GO=int)
1059
1060#endif // HAVE_MUELU_EPETRA
1061
1062
1063
1093 long ExtractNonSerializableData(const Teuchos::ParameterList& inList, Teuchos::ParameterList& serialList, Teuchos::ParameterList& nonSerialList);
1094
1095
1099 void TokenizeStringAndStripWhiteSpace(const std::string & stream, std::vector<std::string> & tokenList, const char* token = ",");
1100
1103 bool IsParamMuemexVariable(const std::string& name);
1104
1107 bool IsParamValidVariable(const std::string& name);
1108
1109#ifdef HAVE_MUELU_EPETRA
1114 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1115 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1116 EpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Epetra_CrsMatrix>& A) {
1117 typedef Xpetra::EpetraCrsMatrixT<GlobalOrdinal, Node> XECrsMatrix;
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 XECrsMatrix(A));
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 EpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Epetra_MultiVector>& V) {
1132 return rcp(new Xpetra::EpetraMultiVectorT<GlobalOrdinal, Node>(V));
1133 }
1134#endif
1135
1136#ifdef HAVE_MUELU_TPETRA
1141 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1142 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1143 TpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >& Atpetra) {
1144 typedef Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> XTCrsMatrix;
1145 typedef Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrix;
1146 typedef Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrixWrap;
1147
1148 RCP<XCrsMatrix> Atmp = rcp(new XTCrsMatrix(Atpetra));
1149 return rcp(new XCrsMatrixWrap(Atmp));
1150 }
1151
1179 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1180 void leftRghtDofScalingWithinNode(const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> & Amat, size_t blkSize, size_t nSweeps, Teuchos::ArrayRCP<Scalar> & rowScaling, Teuchos::ArrayRCP<Scalar> & colScaling) {
1181
1182 LocalOrdinal nBlks = (Amat.getRowMap()->getLocalNumElements())/blkSize;
1183
1184 Teuchos::ArrayRCP<Scalar> rowScaleUpdate(blkSize);
1185 Teuchos::ArrayRCP<Scalar> colScaleUpdate(blkSize);
1186
1187
1188 for (size_t i = 0; i < blkSize; i++) rowScaling[i] = 1.0;
1189 for (size_t i = 0; i < blkSize; i++) colScaling[i] = 1.0;
1190
1191 for (size_t k = 0; k < nSweeps; k++) {
1192 LocalOrdinal row = 0;
1193 for (size_t i = 0; i < blkSize; i++) rowScaleUpdate[i] = 0.0;
1194
1195 for (LocalOrdinal i = 0; i < nBlks; i++) {
1196 for (size_t j = 0; j < blkSize; j++) {
1197 Teuchos::ArrayView<const LocalOrdinal> cols;
1198 Teuchos::ArrayView<const Scalar> vals;
1199 Amat.getLocalRowView(row, cols, vals);
1200
1201 for (size_t kk = 0; kk < Teuchos::as<size_t>(vals.size()); kk++) {
1202 size_t modGuy = (cols[kk]+1)%blkSize;
1203 if (modGuy == 0) modGuy = blkSize;
1204 modGuy--;
1205 rowScaleUpdate[j] += rowScaling[j]*(Teuchos::ScalarTraits<Scalar>::magnitude(vals[kk]))*colScaling[modGuy];
1206 }
1207 row++;
1208 }
1209 }
1210 // combine information across processors
1211 Teuchos::ArrayRCP<Scalar> tempUpdate(blkSize);
1212 Teuchos::reduceAll(*(Amat.getRowMap()->getComm()), Teuchos::REDUCE_SUM, (LocalOrdinal) blkSize, rowScaleUpdate.getRawPtr(), tempUpdate.getRawPtr());
1213 for (size_t i = 0; i < blkSize; i++) rowScaleUpdate[i] = tempUpdate[i];
1214
1215 /* We want to scale by sqrt(1/rowScaleUpdate), but we'll */
1216 /* normalize things by the minimum rowScaleUpdate. That is, the */
1217 /* largest scaling is always one (as normalization is arbitrary).*/
1218
1219 Scalar minUpdate = Teuchos::ScalarTraits<Scalar>::magnitude((rowScaleUpdate[0]/rowScaling[0])/rowScaling[0]);
1220
1221 for (size_t i = 1; i < blkSize; i++) {
1222 Scalar temp = (rowScaleUpdate[i]/rowScaling[i])/rowScaling[i];
1223 if ( Teuchos::ScalarTraits<Scalar>::magnitude(temp) < Teuchos::ScalarTraits<Scalar>::magnitude(minUpdate))
1224 minUpdate = Teuchos::ScalarTraits<Scalar>::magnitude(temp);
1225 }
1226 for (size_t i = 0; i < blkSize; i++) rowScaling[i] *= sqrt(minUpdate / rowScaleUpdate[i]);
1227
1228 row = 0;
1229 for (size_t i = 0; i < blkSize; i++) colScaleUpdate[i] = 0.0;
1230
1231 for (LocalOrdinal i = 0; i < nBlks; i++) {
1232 for (size_t j = 0; j < blkSize; j++) {
1233 Teuchos::ArrayView<const LocalOrdinal> cols;
1234 Teuchos::ArrayView<const Scalar> vals;
1235 Amat.getLocalRowView(row, cols, vals);
1236 for (size_t kk = 0; kk < Teuchos::as<size_t>(vals.size()); kk++) {
1237 size_t modGuy = (cols[kk]+1)%blkSize;
1238 if (modGuy == 0) modGuy = blkSize;
1239 modGuy--;
1240 colScaleUpdate[modGuy] += colScaling[modGuy]* (Teuchos::ScalarTraits<Scalar>::magnitude(vals[kk])) *rowScaling[j];
1241 }
1242 row++;
1243 }
1244 }
1245 Teuchos::reduceAll(*(Amat.getRowMap()->getComm()), Teuchos::REDUCE_SUM, (LocalOrdinal) blkSize, colScaleUpdate.getRawPtr(), tempUpdate.getRawPtr());
1246 for (size_t i = 0; i < blkSize; i++) colScaleUpdate[i] = tempUpdate[i];
1247
1248 /* We want to scale by sqrt(1/colScaleUpdate), but we'll */
1249 /* normalize things by the minimum colScaleUpdate. That is, the */
1250 /* largest scaling is always one (as normalization is arbitrary).*/
1251
1252
1253 minUpdate = Teuchos::ScalarTraits<Scalar>::magnitude((colScaleUpdate[0]/colScaling[0])/colScaling[0]);
1254
1255 for (size_t i = 1; i < blkSize; i++) {
1256 Scalar temp = (colScaleUpdate[i]/colScaling[i])/colScaling[i];
1257 if ( Teuchos::ScalarTraits<Scalar>::magnitude(temp) < Teuchos::ScalarTraits<Scalar>::magnitude(minUpdate))
1258 minUpdate = Teuchos::ScalarTraits<Scalar>::magnitude(temp);
1259 }
1260 for (size_t i = 0; i < blkSize; i++) colScaling[i] *= sqrt(minUpdate/colScaleUpdate[i]);
1261 }
1262 }
1263
1268 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1269 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1270 TpetraFECrs_To_XpetraMatrix(const Teuchos::RCP<Tpetra::FECrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >& Atpetra) {
1271 typedef typename Tpetra::FECrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::crs_matrix_type tpetra_crs_matrix_type;
1272 typedef Xpetra::TpetraCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> XTCrsMatrix;
1273 typedef Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrix;
1274 typedef Xpetra::CrsMatrixWrap<Scalar, LocalOrdinal, GlobalOrdinal, Node> XCrsMatrixWrap;
1275
1276 RCP<XCrsMatrix> Atmp = rcp(new XTCrsMatrix(rcp_dynamic_cast<tpetra_crs_matrix_type>(Atpetra)));
1277 return rcp(new XCrsMatrixWrap(Atmp));
1278 }
1279
1284 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1285 RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1286 TpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >& Vtpetra) {
1287 return rcp(new Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>(Vtpetra));
1288 }
1289
1294 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1295 RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1296 TpetraFEMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Tpetra::FEMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >& Vtpetra) {
1297 typedef Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> MV;
1298 RCP<const MV> Vmv = Teuchos::rcp_dynamic_cast<const MV>(Vtpetra);
1299 return rcp(new Xpetra::TpetraMultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node>(Vmv));
1300 }
1301
1302#endif
1303
1305 template<class T>
1306 std::string toString(const T& what) {
1307 std::ostringstream buf;
1308 buf << what;
1309 return buf.str();
1310 }
1311
1312#ifdef HAVE_MUELU_EPETRA
1317 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1318 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1319 EpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Epetra_CrsMatrix>& A);
1320
1325 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1326 RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1327 EpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Epetra_MultiVector>& V);
1328#endif
1329
1330#ifdef HAVE_MUELU_TPETRA
1335 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1336 RCP<Xpetra::Matrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1337 TpetraCrs_To_XpetraMatrix(const Teuchos::RCP<Tpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> >& Atpetra);
1338
1343 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
1344 RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
1345 TpetraMultiVector_To_XpetraMultiVector(const Teuchos::RCP<Tpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >& Vtpetra);
1346#endif
1347
1348 // Generates a communicator whose only members are other ranks of the baseComm on my node
1349 Teuchos::RCP<const Teuchos::Comm<int> > GenerateNodeComm(RCP<const Teuchos::Comm<int> > & baseComm, int &NodeId, const int reductionFactor);
1350
1351 // Lower case string
1352 std::string lowerCase (const std::string& s);
1353
1354} //namespace MueLu
1355
1356#define MUELU_UTILITIES_SHORT
1357#endif // MUELU_UTILITIES_DECL_HPP
Kokkos::Compat::KokkosSerialWrapperNode EpetraNode
MueLu::DefaultLocalOrdinal LocalOrdinal
MueLu::DefaultScalar Scalar
MueLu::DefaultGlobalOrdinal GlobalOrdinal
MueLu::DefaultNode Node
int NumMyElements() const
Exception indicating invalid cast attempted.
Exception throws to report errors in the internal logical of the program.
Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > Map
static RCP< CrsMatrixWrap > GetThresholdedMatrix(const RCP< Matrix > &Ain, const Scalar threshold, const bool keepDiagonal=true, const GlobalOrdinal expectedNNZperRow=-1)
Threshold a matrix.
static void ZeroDirichletCols(Teuchos::RCP< Matrix > &A, const Teuchos::ArrayRCP< const bool > &dirichletCols, Scalar replaceWith=Teuchos::ScalarTraits< Scalar >::zero())
static RCP< Teuchos::FancyOStream > MakeFancy(std::ostream &os)
static void DetectDirichletColsAndDomains(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Teuchos::ArrayRCP< bool > &dirichletRows, Teuchos::ArrayRCP< bool > dirichletCols, Teuchos::ArrayRCP< bool > dirichletDomain)
Detects Dirichlet columns & domains from a list of Dirichlet rows.
static Teuchos::ScalarTraits< Scalar >::magnitudeType Distance2(const Teuchos::Array< Teuchos::ArrayRCP< const Scalar > > &v, LocalOrdinal i0, LocalOrdinal i1)
Squared distance between two rows in a multivector.
static Teuchos::ArrayRCP< Scalar > GetMatrixDiagonal(const Matrix &A)
Extract Matrix Diagonal.
static void SetRandomSeed(const Teuchos::Comm< int > &comm)
Set seed for random number generator.
static Teuchos::Array< Magnitude > ResidualNorm(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const MultiVector &X, const MultiVector &RHS)
static Teuchos::ArrayRCP< const bool > DetectDirichletRows(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Magnitude &tol=Teuchos::ScalarTraits< Scalar >::zero(), bool count_twos_as_dirichlet=false)
Detect Dirichlet rows.
static Teuchos::ArrayRCP< Magnitude > GetMatrixMaxMinusOffDiagonal(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Return vector containing: max_{i\not=k}(-a_ik), for each for i in the matrix.
static RCP< Vector > GetMatrixDiagonalInverse(const Matrix &A, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps() *100, Scalar valReplacement=Teuchos::ScalarTraits< Scalar >::zero())
Extract Matrix Diagonal.
static Scalar Frobenius(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B)
Frobenius inner product of two matrices.
static Scalar PowerMethod(const Matrix &A, bool scaleByDiag=true, LocalOrdinal niters=10, Magnitude tolerance=1e-2, bool verbose=false, unsigned int seed=123)
Power method.
static void ApplyOAZToMatrixRows(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const std::vector< LocalOrdinal > &dirichletRows)
static RCP< Matrix > Crs2Op(RCP< CrsMatrix > Op)
static Teuchos::ArrayRCP< const bool > DetectDirichletRowsExt(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, bool &bHasZeroDiagonal, const Magnitude &tol=Teuchos::ScalarTraits< Scalar >::zero())
Detect Dirichlet rows (extended version)
static void FindDirichletRows(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, std::vector< LocalOrdinal > &dirichletRows, bool count_twos_as_dirichlet=false)
static void ApplyRowSumCriterion(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Magnitude rowSumTol, Teuchos::ArrayRCP< bool > &dirichletRows)
Apply Rowsum Criterion.
static Teuchos::RCP< Vector > GetLumpedMatrixDiagonal(Matrix const &A, const bool doReciprocal=false, Magnitude tol=Teuchos::ScalarTraits< Scalar >::magnitude(Teuchos::ScalarTraits< Scalar >::zero()), Scalar valReplacement=Teuchos::ScalarTraits< Scalar >::zero(), const bool replaceSingleEntryRowWithZero=false, const bool useAverageAbsDiagVal=false)
Extract Matrix Diagonal of lumped matrix.
static void FindNonZeros(const Teuchos::ArrayRCP< const Scalar > vals, Teuchos::ArrayRCP< bool > nonzeros)
Find non-zero values in an ArrayRCP Compares the value to 2 * machine epsilon.
Xpetra::CrsMatrixWrap< Scalar, LocalOrdinal, GlobalOrdinal, Node > CrsMatrixWrap
static RCP< MultiVector > Residual(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const MultiVector &X, const MultiVector &RHS)
static RCP< Xpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > GetThresholdedGraph(const RCP< Matrix > &A, const Magnitude threshold, const GlobalOrdinal expectedNNZperRow=-1)
Threshold a graph.
static RCP< Vector > GetMatrixOverlappedDiagonal(const Matrix &A)
Extract Overlapped Matrix Diagonal.
static Teuchos::RCP< Vector > GetInverse(Teuchos::RCP< const Vector > v, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps() *100, Scalar valReplacement=Teuchos::ScalarTraits< Scalar >::zero())
Return vector containing inverse of input vector.
static void ZeroDirichletRows(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const std::vector< LocalOrdinal > &dirichletRows, Scalar replaceWith=Teuchos::ScalarTraits< Scalar >::zero())
static const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & MV2TpetraMV(const MultiVector &vec)
static void ApplyRowSumCriterion(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Xpetra::Vector< LocalOrdinal, LocalOrdinal, GlobalOrdinal, Node > &BlockNumber, const Magnitude rowSumTol, Teuchos::ArrayRCP< bool > &dirichletRows)
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)
static RCP< MultiVector > Residual(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &RHS)
static RCP< Vector > GetMatrixOverlappedDiagonal(const Matrix &A)
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)
Xpetra::CrsMatrixWrap< Scalar, LocalOrdinal, GlobalOrdinal, Node > CrsMatrixWrap
static Epetra_MultiVector & MV2NonConstEpetraMV(MultiVector &vec)
static RCP< Vector > GetInverse(Teuchos::RCP< const Vector > v, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps() *100, Scalar tolReplacement=Teuchos::ScalarTraits< Scalar >::zero())
static void SetRandomSeed(const Teuchos::Comm< int > &comm)
Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > Matrix
Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > Map
static RCP< Vector > GetMatrixDiagonalInverse(const Matrix &A, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps() *100, Scalar valReplacement=Teuchos::ScalarTraits< Scalar >::zero())
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 Teuchos::ArrayRCP< const bool > DetectDirichletRows(const Matrix &A, const Magnitude &tol=Teuchos::ScalarTraits< Scalar >::zero(), const bool count_twos_as_dirichlet=false)
static Teuchos::Array< Magnitude > ResidualNorm(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &RHS)
Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > MultiVector
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 Scalar Frobenius(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B)
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)
Xpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > CrsGraph
static const RCP< const Tpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > Map2TpetraMap(const Map &map)
Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > Vector
static RCP< const Epetra_CrsMatrix > Op2EpetraCrs(RCP< const Matrix > Op)
static RCP< CrsMatrixWrap > GetThresholdedMatrix(const RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Scalar threshold, const bool keepDiagonal, const GlobalOrdinal expectedNNZperRow)
static Teuchos::ArrayRCP< const bool > DetectDirichletRowsExt(const Matrix &A, bool &bHasZeroDiagonal, const Magnitude &tol=Teuchos::ScalarTraits< Scalar >::zero())
static const Epetra_MultiVector & MV2EpetraMV(const MultiVector &vec)
static Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > & Op2NonConstTpetraCrs(Matrix &Op)
static Teuchos::ArrayRCP< Magnitude > GetMatrixMaxMinusOffDiagonal(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
static void DetectDirichletColsAndDomains(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Teuchos::ArrayRCP< bool > &dirichletRows, Teuchos::ArrayRCP< bool > dirichletCols, Teuchos::ArrayRCP< bool > dirichletDomains)
static void ApplyRowSumCriterion(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Magnitude rowSumTol, Teuchos::ArrayRCP< bool > &dirichletRows)
static RCP< const Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraCrs(RCP< const Matrix > Op)
Xpetra::EpetraMultiVectorT< GlobalOrdinal, Node > EpetraMultiVector
static Teuchos::Array< Magnitude > ResidualNorm(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &RHS, Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Resid)
static Scalar PowerMethod(const Matrix &A, const RCP< Vector > &invDiag, LocalOrdinal niters=10, Magnitude tolerance=1e-2, bool verbose=false, unsigned int seed=123)
static Teuchos::RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > GetLumpedMatrixDiagonal(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > const &A, const bool doReciprocal=false, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps() *100, Scalar tolReplacement=Teuchos::ScalarTraits< Scalar >::zero(), const bool replaceSingleEntryRowWithZero=false, const bool useAverageAbsDiagVal=false)
static RCP< Teuchos::FancyOStream > MakeFancy(std::ostream &os)
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< CrsGraph > GetThresholdedGraph(const RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Magnitude threshold, const GlobalOrdinal expectedNNZperRow)
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 Teuchos::ScalarTraits< Scalar >::magnitudeType Distance2(const Teuchos::Array< Teuchos::ArrayRCP< const Scalar > > &v, LocalOrdinal i0, LocalOrdinal i1)
static RCP< Epetra_CrsMatrix > Op2NonConstEpetraCrs(RCP< Matrix > Op)
static void Residual(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &RHS, Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Resid)
static Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > & MV2NonConstTpetraMV(MultiVector &vec)
static const Epetra_CrsMatrix & Op2EpetraCrs(const Matrix &Op)
Xpetra::TpetraOperator< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraOperator
static RCP< Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraBlockCrs(RCP< Matrix > Op)
static Teuchos::ScalarTraits< Scalar >::magnitudeType Distance2(const Teuchos::ArrayView< double > &weight, const Teuchos::Array< Teuchos::ArrayRCP< const Scalar > > &v, LocalOrdinal i0, LocalOrdinal i1)
Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > CrsMatrix
Xpetra::EpetraCrsMatrixT< GlobalOrdinal, Node > EpetraCrsMatrix
static void FindNonZeros(const Teuchos::ArrayRCP< const Scalar > vals, Teuchos::ArrayRCP< bool > nonzeros)
static Teuchos::ArrayRCP< Magnitude > GetMatrixMaxMinusOffDiagonal(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Xpetra::Vector< LocalOrdinal, LocalOrdinal, GlobalOrdinal, Node > &BlockNumber)
Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > Operator
static Teuchos::ArrayRCP< Scalar > GetMatrixDiagonal(const Matrix &A)
static Scalar PowerMethod(const Matrix &A, bool scaleByDiag=true, LocalOrdinal niters=10, Magnitude tolerance=1e-2, bool verbose=false, unsigned int seed=123)
MueLu utility class.
static Teuchos::ArrayRCP< Scalar > GetMatrixDiagonal(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
static void ZeroDirichletRows(Teuchos::RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &X, const Teuchos::ArrayRCP< const bool > &dirichletRows, Scalar replaceWith=Teuchos::ScalarTraits< Scalar >::zero())
static Teuchos::ArrayRCP< const bool > DetectDirichletRowsExt(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, bool &bHasZeroDiagonal, const Magnitude &tol=Teuchos::ScalarTraits< Scalar >::zero())
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)
static RCP< Teuchos::FancyOStream > MakeFancy(std::ostream &os)
Teuchos::ScalarTraits< Scalar >::magnitudeType Magnitude
static void MyOldScaleMatrix_Tpetra(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
static void ZeroDirichletCols(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Teuchos::ArrayRCP< const bool > &dirichletCols, Scalar replaceWith=Teuchos::ScalarTraits< Scalar >::zero())
static Teuchos::Array< Magnitude > ResidualNorm(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &RHS, Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Resid)
static void Residual(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &RHS, Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Resid)
static void FindNonZeros(const Teuchos::ArrayRCP< const Scalar > vals, Teuchos::ArrayRCP< bool > nonzeros)
static Teuchos::RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > GetLumpedMatrixDiagonal(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > const &A, const bool doReciprocal=false, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps() *100, Scalar tolReplacement=Teuchos::ScalarTraits< Scalar >::zero(), const bool replaceSingleEntryRowWithZero=false, const bool useAverageAbsDiagVal=false)
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< Xpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > GetThresholdedGraph(const RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Magnitude threshold, const GlobalOrdinal expectedNNZperRow)
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 Teuchos::Array< Magnitude > ResidualNorm(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &RHS)
static Scalar PowerMethod(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, bool scaleByDiag=true, LocalOrdinal niters=10, Magnitude tolerance=1e-2, bool verbose=false, unsigned int seed=123)
static Teuchos::ArrayRCP< Magnitude > GetMatrixMaxMinusOffDiagonal(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Xpetra::Vector< LocalOrdinal, LocalOrdinal, GlobalOrdinal, Node > &BlockNumber)
static const Epetra_Map & Map2EpetraMap(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &map)
static void ApplyRowSumCriterion(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Magnitude rowSumTol, Teuchos::ArrayRCP< bool > &dirichletRows)
static Teuchos::ArrayRCP< Magnitude > GetMatrixMaxMinusOffDiagonal(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
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 void FindDirichletRows(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, std::vector< LocalOrdinal > &dirichletRows, bool count_twos_as_dirichlet=false)
static void ZeroDirichletRows(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const std::vector< LocalOrdinal > &dirichletRows, Scalar replaceWith=Teuchos::ScalarTraits< Scalar >::zero())
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 Teuchos::ScalarTraits< Scalar >::magnitudeType Distance2(const Teuchos::ArrayView< double > &weight, const Teuchos::Array< Teuchos::ArrayRCP< const Scalar > > &v, LocalOrdinal i0, LocalOrdinal i1)
static void MyOldScaleMatrix_Epetra(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< Scalar > &scalingVector, bool doFillComplete, bool doOptimizeStorage)
static RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > GetMatrixDiagonalInverse(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps() *100, Scalar valReplacement=Teuchos::ScalarTraits< Scalar >::zero())
static Scalar Frobenius(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &B)
static void ApplyOAZToMatrixRows(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const std::vector< LocalOrdinal > &dirichletRows)
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< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Residual(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &RHS)
static void ZeroDirichletRows(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Teuchos::ArrayRCP< const bool > &dirichletRows, Scalar replaceWith=Teuchos::ScalarTraits< Scalar >::zero())
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< Xpetra::CrsMatrixWrap< Scalar, LocalOrdinal, GlobalOrdinal, Node > > GetThresholdedMatrix(const RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Scalar threshold, const bool keepDiagonal, const GlobalOrdinal expectedNNZperRow)
static RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > GetMatrixOverlappedDiagonal(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
static RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraBlockCrs(RCP< const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static Scalar PowerMethod(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &invDiag, LocalOrdinal niters=10, Magnitude tolerance=1e-2, bool verbose=false, unsigned int seed=123)
static RCP< const Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2TpetraRow(RCP< const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static void ApplyOAZToMatrixRows(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Teuchos::ArrayRCP< const bool > &dirichletRows)
static RCP< Xpetra::MultiVector< typename Teuchos::ScalarTraits< Scalar >::magnitudeType, LocalOrdinal, GlobalOrdinal, Node > > ExtractCoordinatesFromParameterList(ParameterList &paramList)
static void ApplyRowSumCriterion(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Xpetra::Vector< LocalOrdinal, LocalOrdinal, GlobalOrdinal, Node > &BlockNumber, const Magnitude rowSumTol, Teuchos::ArrayRCP< bool > &dirichletRows)
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 void SetRandomSeed(const Teuchos::Comm< int > &comm)
static Teuchos::ScalarTraits< Scalar >::magnitudeType Distance2(const Teuchos::Array< Teuchos::ArrayRCP< const Scalar > > &v, LocalOrdinal i0, LocalOrdinal i1)
static Teuchos::RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > GetInverse(Teuchos::RCP< const Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > v, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps() *100, Scalar tolReplacement=Teuchos::ScalarTraits< Scalar >::zero())
static RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Crs2Op(RCP< Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static void DetectDirichletColsAndDomains(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Teuchos::ArrayRCP< bool > &dirichletRows, Teuchos::ArrayRCP< bool > dirichletCols, Teuchos::ArrayRCP< bool > dirichletDomains)
static Teuchos::ArrayRCP< const bool > DetectDirichletRows(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Magnitude &tol=Teuchos::ScalarTraits< Scalar >::magnitude(0.), const bool count_twos_as_dirichlet=false)
Namespace for MueLu classes and methods.
bool IsParamMuemexVariable(const std::string &name)
KokkosClassic::DefaultNode::DefaultNodeType DefaultNode
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)
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)
magnitude_type tolerance