46#ifndef MUELU_UTILITIESBASE_DECL_HPP
47#define MUELU_UTILITIESBASE_DECL_HPP
57#include <Teuchos_DefaultComm.hpp>
58#include <Teuchos_ScalarTraits.hpp>
59#include <Teuchos_ParameterList.hpp>
61#include <Xpetra_BlockedCrsMatrix_fwd.hpp>
62#include <Xpetra_CrsGraphFactory_fwd.hpp>
63#include <Xpetra_CrsGraph_fwd.hpp>
64#include <Xpetra_CrsMatrix_fwd.hpp>
65#include <Xpetra_CrsMatrixWrap_fwd.hpp>
66#include <Xpetra_Map_fwd.hpp>
67#include <Xpetra_BlockedMap_fwd.hpp>
68#include <Xpetra_MapFactory_fwd.hpp>
69#include <Xpetra_Matrix_fwd.hpp>
70#include <Xpetra_MatrixFactory_fwd.hpp>
71#include <Xpetra_MultiVector_fwd.hpp>
72#include <Xpetra_MultiVectorFactory_fwd.hpp>
73#include <Xpetra_Operator_fwd.hpp>
74#include <Xpetra_Vector_fwd.hpp>
75#include <Xpetra_BlockedMultiVector.hpp>
76#include <Xpetra_BlockedVector.hpp>
77#include <Xpetra_VectorFactory_fwd.hpp>
78#include <Xpetra_ExportFactory.hpp>
80#include <Xpetra_Import.hpp>
81#include <Xpetra_ImportFactory.hpp>
82#include <Xpetra_MatrixMatrix.hpp>
83#include <Xpetra_CrsGraph.hpp>
84#include <Xpetra_CrsGraphFactory.hpp>
85#include <Xpetra_CrsMatrixWrap.hpp>
86#include <Xpetra_StridedMap.hpp>
94#define MueLu_sumAll(rcpComm, in, out) \
95 Teuchos::reduceAll(*rcpComm, Teuchos::REDUCE_SUM, in, Teuchos::outArg(out))
96#define MueLu_minAll(rcpComm, in, out) \
97 Teuchos::reduceAll(*rcpComm, Teuchos::REDUCE_MIN, in, Teuchos::outArg(out))
98#define MueLu_maxAll(rcpComm, in, out) \
99 Teuchos::reduceAll(*rcpComm, Teuchos::REDUCE_MAX, in, Teuchos::outArg(out))
114#undef MUELU_UTILITIESBASE_SHORT
117 using CrsGraph = Xpetra::CrsGraph<LocalOrdinal,GlobalOrdinal,Node>;
119 using CrsMatrixWrap = Xpetra::CrsMatrixWrap<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
120 using CrsMatrix = Xpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
121 using Matrix = Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
122 using Vector = Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
123 using BlockedVector = Xpetra::BlockedVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
124 using MultiVector = Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>;
126 using BlockedMap = Xpetra::BlockedMap<LocalOrdinal,GlobalOrdinal,Node>;
127 using Map = Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node>;
129 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType
Magnitude;
134 return Teuchos::null;
152 size_t nnz =
Ain->getNumEntriesInLocalRow(
row);
154 Teuchos::ArrayView<const LocalOrdinal> indices;
155 Teuchos::ArrayView<const Scalar>
vals;
160 Teuchos::ArrayRCP<GlobalOrdinal>
indout(indices.size(),Teuchos::ScalarTraits<GlobalOrdinal>::zero());
161 Teuchos::ArrayRCP<Scalar>
valout(indices.size(),Teuchos::ScalarTraits<Scalar>::zero());
166 for(
size_t i=0; i<(
size_t)indices.size(); i++) {
167 if(Teuchos::ScalarTraits<Scalar>::magnitude(
vals[i]) > Teuchos::ScalarTraits<Scalar>::magnitude(
threshold) || indices[i]==
lclColIdx) {
174 for(
size_t i=0; i<(
size_t)indices.size(); i++) {
175 if(Teuchos::ScalarTraits<Scalar>::magnitude(
vals[i]) > Teuchos::ScalarTraits<Scalar>::magnitude(
threshold)) {
187 Aout->fillComplete(
Ain->getDomainMap(),
Ain->getRangeMap());
200 using STS = Teuchos::ScalarTraits<Scalar>;
206 for(
size_t row=0;
row<A->getLocalNumRows();
row++)
210 A->getLocalRowView(
row, indices,
vals);
215 const Scalar Dk = STS::magnitude(D[col]) > 0.0 ? STS::magnitude(D[col]) : 1.0;
218 for(
size_t i=0; i<
size_t(indices.size()); i++)
220 if(col == indices[i] || STS::magnitude(STS::squareroot(
Dk)*
vals[i]*STS::squareroot(
Dk)) > STS::magnitude(
threshold))
221 indicesNew.append(A->getColMap()->getGlobalElement(indices[i]));
237 size_t numRows = A.getRowMap()->getLocalNumElements();
238 Teuchos::ArrayRCP<Scalar> diag(numRows);
239 Teuchos::ArrayView<const LocalOrdinal>
cols;
240 Teuchos::ArrayView<const Scalar>
vals;
241 for (
size_t i = 0; i < numRows; ++i) {
244 for (; j <
cols.size(); ++j) {
245 if (Teuchos::as<size_t>(
cols[j]) == i) {
250 if (j ==
cols.size()) {
252 diag[i] = Teuchos::ScalarTraits<Scalar>::zero();
265 Teuchos::TimeMonitor
MM = *Teuchos::TimeMonitor::getNewTimer(
"UtilitiesBase::GetMatrixDiagonalInverse");
268 RCP<Vector> diag = Xpetra::VectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(rowMap,
true);
270 A.getLocalDiagCopy(*diag);
284 Magnitude tol = Teuchos::ScalarTraits<Scalar>::magnitude(Teuchos::ScalarTraits<Scalar>::zero()),
289 typedef Teuchos::ScalarTraits<Scalar> TST;
292 const Scalar zero = TST::zero();
293 const Scalar one = TST::one();
296 Teuchos::RCP<const Matrix>
rcpA = Teuchos::rcpFromRef(A);
299 Teuchos::rcp_dynamic_cast<const Xpetra::BlockedCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(
rcpA);
300 if(
bA == Teuchos::null) {
302 diag = Xpetra::VectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(rowMap,
true);
304 Teuchos::Array<Scalar>
regSum(diag->getLocalLength());
305 Teuchos::ArrayView<const LocalOrdinal>
cols;
306 Teuchos::ArrayView<const Scalar>
vals;
308 std::vector<int>
nnzPerRow(rowMap->getLocalNumElements());
316 for (
size_t i = 0; i < rowMap->getLocalNumElements(); ++i) {
326 if (
static_cast<size_t>(
cols[j]) == i)
335 for (
size_t i = 0; i < rowMap->getLocalNumElements(); ++i) {
351 "UtilitiesBase::GetLumpedMatrixDiagonal(): extracting reciprocal of diagonal of a blocked matrix is not supported");
352 diag = Xpetra::VectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(
bA->getRangeMapExtractor()->getFullMap(),
true);
355 for (
size_t col = 0; col <
bA->Cols(); ++col) {
356 if (!
bA->getMatrix(
row,col).is_null()) {
358 bool bThyraMode =
bA->getRangeMapExtractor()->getThyraMode() && (Teuchos::rcp_dynamic_cast<Xpetra::BlockedCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(
bA->getMatrix(
row,col)) == Teuchos::null);
361 ddtemp->update(Teuchos::as<Scalar>(1.0),*
dd,Teuchos::as<Scalar>(1.0));
379 size_t numRows = A.getRowMap()->getLocalNumElements();
381 Teuchos::ArrayRCP<Magnitude>
maxvec(numRows);
382 Teuchos::ArrayView<const LocalOrdinal>
cols;
383 Teuchos::ArrayView<const Scalar>
vals;
384 for (
size_t i = 0; i < numRows; ++i) {
388 if (Teuchos::as<size_t>(
cols[j]) != i) {
389 mymax = std::max(
mymax,-Teuchos::ScalarTraits<Scalar>::real(
vals[j]));
397 static Teuchos::ArrayRCP<Magnitude>
GetMatrixMaxMinusOffDiagonal(
const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A,
const Xpetra::Vector<LocalOrdinal,LocalOrdinal,GlobalOrdinal,Node> &
BlockNumber) {
402 size_t numRows = A.getRowMap()->getLocalNumElements();
404 Teuchos::ArrayRCP<Magnitude>
maxvec(numRows);
405 Teuchos::ArrayView<const LocalOrdinal>
cols;
406 Teuchos::ArrayView<const Scalar>
vals;
407 for (
size_t i = 0; i < numRows; ++i) {
412 mymax = std::max(
mymax,-Teuchos::ScalarTraits<Scalar>::real(
vals[j]));
431 RCP<Vector> ret = Xpetra::VectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(v->getMap(),
true);
435 if(
bv.is_null() ==
false) {
439 for(
size_t r = 0;
r <
bmap->getNumMaps(); ++
r) {
452 if(Teuchos::ScalarTraits<Scalar>::magnitude(
inputVals[i]) >
tol)
474 RCP<Vector> localDiag = Xpetra::VectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(rowMap);
480 Teuchos::ArrayRCP<size_t> offsets;
481 crsOp->getLocalDiagOffsets(offsets);
492 RCP<Vector> diagonal = Xpetra::VectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(colMap);
494 importer = A.getCrsGraph()->getImporter();
496 importer = Xpetra::ImportFactory<LocalOrdinal,GlobalOrdinal,Node>::Build(rowMap, colMap);
514 using STS =
typename Teuchos::ScalarTraits<SC>;
526 size_t nnz = A.getNumEntriesInLocalRow(
row);
529 A.getLocalRowView(
row, indices,
vals);
542 importer = A.getCrsGraph()->getImporter();
544 importer = Xpetra::ImportFactory<LO,GO,Node>::Build(rowMap, colMap);
552 static RCP<Xpetra::Vector<Magnitude,LocalOrdinal,GlobalOrdinal,Node> >
555 using STS =
typename Teuchos::ScalarTraits<Scalar>;
556 using MTS =
typename Teuchos::ScalarTraits<Magnitude>;
572 size_t nnz = A.getNumEntriesInLocalRow(
rowIdx);
588 importer = A.getCrsGraph()->getImporter();
590 importer = Xpetra::ImportFactory<LO,GO,Node>::Build(rowMap, colMap);
604 const size_t numVecs = X.getNumVectors();
613 const size_t numVecs = X.getNumVectors();
622 const size_t numVecs = X.getNumVectors();
624 RCP<MultiVector> RES = Xpetra::MultiVectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(RHS.getMap(),
numVecs,
false);
625 Op.residual(X,RHS,*
RES);
633 Op.residual(X,RHS,
Resid);
651 "Utils::PowerMethod: operator must have domain and range maps that are equivalent.");
656 RCP<Vector> diagVec = Xpetra::VectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(A.getRowMap());
658 diagInvVec = Xpetra::VectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(A.getRowMap());
680 "Utils::PowerMethod: operator must have domain and range maps that are equivalent.");
683 RCP<Vector> q = Xpetra::VectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(A.getDomainMap());
684 RCP<Vector> r = Xpetra::VectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(A.getRangeMap());
685 RCP<Vector> z = Xpetra::VectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(A.getRangeMap());
690 Teuchos::Array<Magnitude>
norms(1);
692 typedef Teuchos::ScalarTraits<Scalar> STS;
694 const Scalar zero = STS::zero(), one = STS::one();
697 Magnitude residual = STS::magnitude(zero);
702 q->update(one/
norms[0], *z, zero);
705 z->elementWiseMultiply(one, *
diagInvVec, *z, zero);
709 r->update(1.0, *z, -
lambda, *
q, zero);
713 std::cout <<
"Iter = " <<
iter
715 <<
" Residual of A*q - lambda*q = " << residual
737 Scalar d = Teuchos::ScalarTraits<Scalar>::zero();
739 d += (v[j][
i0] - v[j][
i1])*(v[j][
i0] - v[j][
i1]);
741 return Teuchos::ScalarTraits<Scalar>::magnitude(
d);
750 using MT =
typename Teuchos::ScalarTraits<Scalar>::magnitudeType;
752 Scalar d = Teuchos::ScalarTraits<Scalar>::zero();
754 d += Teuchos::as<MT>(
weight[j])*(v[j][
i0] - v[j][
i1])*(v[j][
i0] - v[j][
i1]);
756 return Teuchos::ScalarTraits<Scalar>::magnitude(
d);
774 typedef Teuchos::ScalarTraits<Scalar> STS;
780 A.getLocalRowView(
row, indices,
vals);
781 size_t nnz = A.getNumEntriesInLocalRow(
row);
784 for (col = 0; col < nnz; col++)
785 if ( (indices[col] !=
row) && STS::magnitude(
vals[col]) >
tol) {
798 A.getLocalRowView(
row, indices,
vals);
799 size_t nnz = A.getNumEntriesInLocalRow(
row);
801 for (
size_t col = 0; col < nnz; col++)
802 if ( (indices[col] !=
row) && STS::magnitude(
vals[col]) >
tol) {
829 Teuchos::RCP<Vector>
diagVec = Xpetra::VectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(A.getRowMap());
834 typedef Teuchos::ScalarTraits<Scalar> STS;
839 A.getLocalRowView(
row, indices,
vals);
842 for (
decltype(indices.size()) col = 0; col < indices.size(); col++) {
843 if ( indices[col] !=
row) {
865 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType magnitudeType;
866 const magnitudeType eps = 2.0*Teuchos::ScalarTraits<magnitudeType>::eps();
868 nonzeros[i] = (Teuchos::ScalarTraits<Scalar>::magnitude(
vals[i]) > eps);
884 const Scalar one = Teuchos::ScalarTraits<Scalar>::one();
897 A.getLocalRowView(i,indices,values);
906 globalColsToZero = Xpetra::MultiVectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(
domMap, 1,
true);
933 typedef Teuchos::ScalarTraits<Scalar> STS;
934 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType MT;
935 typedef Teuchos::ScalarTraits<MT>
MTS;
938 size_t nnz = A.getNumEntriesInLocalRow(
row);
941 A.getLocalRowView(
row, indices,
vals);
961 typedef Teuchos::ScalarTraits<Scalar> STS;
962 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType MT;
963 typedef Teuchos::ScalarTraits<MT>
MTS;
970 size_t nnz = A.getNumEntriesInLocalRow(
row);
973 A.getLocalRowView(
row, indices,
vals);
1005 static Teuchos::ArrayRCP<const bool>
DetectDirichletCols(
const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A,
1007 Scalar zero = Teuchos::ScalarTraits<Scalar>::zero();
1008 Scalar one = Teuchos::ScalarTraits<Scalar>::one();
1009 Teuchos::RCP<const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> >
domMap = A.getDomainMap();
1010 Teuchos::RCP<const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > colMap = A.getColMap();
1011 Teuchos::RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
myColsToZero = Xpetra::MultiVectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(colMap,1);
1016 Teuchos::ArrayView<const LocalOrdinal> indices;
1017 Teuchos::ArrayView<const Scalar> values;
1018 A.getLocalRowView(i,indices,values);
1024 Teuchos::RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
globalColsToZero = Xpetra::MultiVectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(
domMap,1);
1026 Teuchos::RCP<Xpetra::Export<LocalOrdinal,GlobalOrdinal,Node> >
exporter = Xpetra::ExportFactory<LocalOrdinal,GlobalOrdinal,Node>::Build(colMap,
domMap);
1032 Teuchos::ArrayRCP<bool>
dirichletCols(colMap->getLocalNumElements(),
true);
1033 Magnitude eps = Teuchos::ScalarTraits<Magnitude>::eps();
1034 for(
size_t i=0; i<colMap->getLocalNumElements(); i++) {
1045 static Scalar Frobenius(
const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A,
const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& B) {
1056 Teuchos::ArrayView<const LocalOrdinal>
indA,
indB;
1057 Teuchos::ArrayView<const Scalar>
valA,
valB;
1073 Scalar zero = Teuchos::ScalarTraits<Scalar> ::zero(),
f = zero,
gf;
1074 size_t numRows = A.getLocalNumRows();
1075 for (
size_t i = 0; i < numRows; i++) {
1082 for (
size_t j = 0; j <
nnzB; j++)
1085 for (
size_t j = 0; j <
nnzA; j++) {
1094 for (
size_t j = 0; j <
nnzB; j++)
1121 std::ostringstream
errStr;
1122 errStr <<
"Error detected with random seed = " <<
mySeed <<
". It should be in the interval [1,2^31-2].";
1127 Teuchos::ScalarTraits<Scalar>::seedrandom(
mySeed);
1139 static void FindDirichletRows(Teuchos::RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> > &A,
1141 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType MT;
1143 for(
size_t i=0; i<A->getLocalNumRows(); i++) {
1144 Teuchos::ArrayView<const LocalOrdinal> indices;
1145 Teuchos::ArrayView<const Scalar> values;
1146 A->getLocalRowView(i,indices,values);
1148 for (
size_t j=0; j<(
size_t)indices.size(); j++) {
1149 if (Teuchos::ScalarTraits<Scalar>::magnitude(values[j]) > Teuchos::ScalarTraits<MT>::eps()) {
1165 Scalar one =Teuchos::ScalarTraits<Scalar>::one();
1166 Scalar zero =Teuchos::ScalarTraits<Scalar>::zero();
1171 Teuchos::ArrayView<const LocalOrdinal> indices;
1172 Teuchos::ArrayView<const Scalar> values;
1176 for(
size_t j=0; j<(
size_t)indices.size(); j++) {
1195 const Scalar one = Teuchos::ScalarTraits<Scalar>::one();
1196 const Scalar zero = Teuchos::ScalarTraits<Scalar>::zero();
1202 Teuchos::ArrayView<const LocalOrdinal> indices;
1203 Teuchos::ArrayView<const Scalar> values;
1204 A->getLocalRowView(i,indices,values);
1206 Teuchos::ArrayRCP<Scalar>
valuesNC(values.size());
1207 for(
size_t j=0; j<(
size_t)indices.size(); j++) {
1213 A->replaceLocalValues(i,indices,
valuesNC());
1221 static void ZeroDirichletRows(Teuchos::RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& A,
1225 Teuchos::ArrayView<const LocalOrdinal> indices;
1226 Teuchos::ArrayView<const Scalar> values;
1230 for(
size_t j=0; j<(
size_t)indices.size(); j++)
1237 static void ZeroDirichletRows(Teuchos::RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& A,
1243 Teuchos::ArrayView<const LocalOrdinal> indices;
1244 Teuchos::ArrayView<const Scalar> values;
1245 A->getLocalRowView(i,indices,values);
1248 for(
size_t j=0; j<(
size_t)indices.size(); j++)
1256 static void ZeroDirichletRows(Teuchos::RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& X,
1262 for(
size_t j=0; j<X->getNumVectors(); j++)
1274 for(
size_t i=0; i<A->getLocalNumRows(); i++) {
1275 Teuchos::ArrayView<const LocalOrdinal> indices;
1276 Teuchos::ArrayView<const Scalar> values;
1277 A->getLocalRowView(i,indices,values);
1288 Teuchos::RCP<Xpetra::Vector<int,LocalOrdinal,GlobalOrdinal,Node> >&
isDirichletRow,
1289 Teuchos::RCP<Xpetra::Vector<int,LocalOrdinal,GlobalOrdinal,Node> >&
isDirichletCol) {
1292 if(!A->getRowMap()->isSameAs(*A->getDomainMap())) {
1293 throw std::runtime_error(
"UtilitiesBase::FindDirichletRowsAndPropagateToCols row and domain maps must match.");
1303 printf(
"[%d] DirichletRow Ids = ",A->getRowMap()->getComm()->getRank());
1310 isDirichletRow = Xpetra::VectorFactory<int,LocalOrdinal,GlobalOrdinal,Node>::Build(A->getRowMap(),
true);
1311 isDirichletCol = Xpetra::VectorFactory<int,LocalOrdinal,GlobalOrdinal,Node>::Build(A->getColMap(),
true);
1315 Teuchos::ArrayView<int>
dr =
dr_rcp();
1317 Teuchos::ArrayView<int>
dc =
dc_rcp();
1333 const Xpetra::Import<LocalOrdinal,GlobalOrdinal,Node> & Importer) {
1334 typedef Xpetra::Vector<int,LocalOrdinal,GlobalOrdinal,Node>
IntVector;
1344 if(!Importer.getSourceMap()->isCompatible(*
fullMap))
1345 throw std::runtime_error(
"GenerateBlockedTargetMap(): Map compatibility error");
1353 for(
size_t j=0; j<
map->getLocalNumElements(); j++) {
1364 Teuchos::ArrayView<const int> data =
dataRCP();
1369 for(
size_t i=0; i<
targetMap->getLocalNumElements(); i++) {
1386 static bool MapsAreNested(
const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node>& rowMap,
const Xpetra::Map<LocalOrdinal,GlobalOrdinal,Node>& colMap) {
1415#define MUELU_UTILITIESBASE_SHORT
#define MueLu_sumAll(rcpComm, in, out)
MueLu::DefaultLocalOrdinal LocalOrdinal
MueLu::DefaultScalar Scalar
MueLu::DefaultGlobalOrdinal GlobalOrdinal
Exception throws to report incompatible objects (like maps).
Exception throws to report errors in the internal logical of the program.
Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > MultiVector
static RCP< Vector > GetMatrixOverlappedDeletedRowsum(const Matrix &A)
Extract Overlapped Matrix Deleted Rowsum.
Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > Map
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 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.
Teuchos::ScalarTraits< Scalar >::magnitudeType Magnitude
Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > CrsMatrix
static Scalar PowerMethod(const Matrix &A, const RCP< Vector > &diagInvVec, LocalOrdinal niters=10, Magnitude tolerance=1e-2, bool verbose=false, unsigned int seed=123)
Power method.
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 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)
Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > Vector
static Teuchos::Array< Magnitude > ResidualNorm(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const MultiVector &X, const MultiVector &RHS)
Xpetra::BlockedVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > BlockedVector
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 RCP< const Xpetra::BlockedMap< LocalOrdinal, GlobalOrdinal, Node > > GeneratedBlockedTargetMap(const Xpetra::BlockedMap< LocalOrdinal, GlobalOrdinal, Node > &sourceBlockedMap, const Xpetra::Import< LocalOrdinal, GlobalOrdinal, Node > &Importer)
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 ZeroDirichletRows(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Teuchos::ArrayRCP< const bool > &dirichletRows, Scalar replaceWith=Teuchos::ScalarTraits< Scalar >::zero())
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.
static bool MapsAreNested(const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &rowMap, const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &colMap)
Xpetra::CrsGraphFactory< LocalOrdinal, GlobalOrdinal, Node > CrsGraphFactory
static void Residual(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const MultiVector &X, const MultiVector &RHS, MultiVector &Resid)
static Teuchos::ScalarTraits< Scalar >::magnitudeType Distance2(const Teuchos::ArrayView< double > &weight, const Teuchos::Array< Teuchos::ArrayRCP< const Scalar > > &v, LocalOrdinal i0, LocalOrdinal i1)
Weighted squared distance between two rows in a multivector.
Xpetra::CrsMatrixWrap< Scalar, LocalOrdinal, GlobalOrdinal, Node > CrsMatrixWrap
static Teuchos::ArrayRCP< const bool > DetectDirichletCols(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Teuchos::ArrayRCP< const bool > &dirichletRows)
Detect Dirichlet columns based on Dirichlet rows.
static RCP< Xpetra::Vector< Magnitude, LocalOrdinal, GlobalOrdinal, Node > > GetMatrixOverlappedAbsDeletedRowsum(const Matrix &A)
static Teuchos::Array< Magnitude > ResidualNorm(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const MultiVector &X, const MultiVector &RHS, MultiVector &Resid)
static Teuchos::ArrayRCP< Magnitude > GetMatrixMaxMinusOffDiagonal(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Xpetra::Vector< LocalOrdinal, LocalOrdinal, GlobalOrdinal, Node > &BlockNumber)
Xpetra::BlockedMultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > BlockedMultiVector
Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > Matrix
Xpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > CrsGraph
static RCP< MultiVector > Residual(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const MultiVector &X, const MultiVector &RHS)
Xpetra::BlockedMap< LocalOrdinal, GlobalOrdinal, Node > BlockedMap
static RCP< Xpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > GetThresholdedGraph(const RCP< Matrix > &A, const Magnitude threshold, const GlobalOrdinal expectedNNZperRow=-1)
Threshold a graph.
static void ApplyOAZToMatrixRows(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Teuchos::ArrayRCP< const bool > &dirichletRows)
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 FindDirichletRowsAndPropagateToCols(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, Teuchos::RCP< Xpetra::Vector< int, LocalOrdinal, GlobalOrdinal, Node > > &isDirichletRow, Teuchos::RCP< Xpetra::Vector< int, LocalOrdinal, GlobalOrdinal, Node > > &isDirichletCol)
static void ZeroDirichletRows(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const std::vector< LocalOrdinal > &dirichletRows, Scalar replaceWith=Teuchos::ScalarTraits< Scalar >::zero())
Namespace for MueLu classes and methods.
KokkosClassic::DefaultNode::DefaultNodeType DefaultNode