46 #ifndef MUELU_UTILITIES_KOKKOS_DEF_HPP 47 #define MUELU_UTILITIES_KOKKOS_DEF_HPP 49 #include <Teuchos_DefaultComm.hpp> 50 #include <Teuchos_ParameterList.hpp> 54 #ifdef HAVE_MUELU_EPETRA 56 # include "Epetra_MpiComm.h" 60 #include <Kokkos_ArithTraits.hpp> 61 #include <Kokkos_Core.hpp> 62 #include <KokkosSparse_CrsMatrix.hpp> 64 #if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_EPETRAEXT) 71 #include <Xpetra_EpetraUtils.hpp> 72 #include <Xpetra_EpetraMultiVector.hpp> 76 #ifdef HAVE_MUELU_TPETRA 77 #include <MatrixMarket_Tpetra.hpp> 78 #include <Tpetra_RowMatrixTransposer.hpp> 79 #include <TpetraExt_MatrixMatrix.hpp> 80 #include <Xpetra_TpetraMultiVector.hpp> 81 #include <Xpetra_TpetraCrsMatrix.hpp> 82 #include <Xpetra_TpetraBlockCrsMatrix.hpp> 85 #ifdef HAVE_MUELU_EPETRA 86 #include <Xpetra_EpetraMap.hpp> 89 #include <Xpetra_BlockedCrsMatrix.hpp> 90 #include <Xpetra_DefaultPlatform.hpp> 91 #include <Xpetra_Import.hpp> 92 #include <Xpetra_ImportFactory.hpp> 93 #include <Xpetra_Map.hpp> 94 #include <Xpetra_MapFactory.hpp> 95 #include <Xpetra_Matrix.hpp> 96 #include <Xpetra_MatrixMatrix.hpp> 97 #include <Xpetra_MatrixFactory.hpp> 98 #include <Xpetra_MultiVector.hpp> 99 #include <Xpetra_MultiVectorFactory.hpp> 100 #include <Xpetra_Operator.hpp> 101 #include <Xpetra_Vector.hpp> 102 #include <Xpetra_VectorFactory.hpp> 108 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
109 Teuchos::ArrayRCP<Scalar> Utilities_kokkos<Scalar, LocalOrdinal, GlobalOrdinal, Node>::GetMatrixDiagonal(
const Matrix& A) {
112 size_t numRows = A.getRowMap()->getNodeNumElements();
113 Teuchos::ArrayRCP<SC> diag(numRows);
115 Teuchos::ArrayView<const LO> cols;
116 Teuchos::ArrayView<const SC> vals;
117 for (
size_t i = 0; i < numRows; ++i) {
118 A.getLocalRowView(i, cols, vals);
121 for (; j < cols.size(); ++j) {
122 if (Teuchos::as<size_t>(cols[j]) == i) {
127 if (j == cols.size()) {
129 diag[i] = Teuchos::ScalarTraits<SC>::zero();
136 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
137 Teuchos::RCP<Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Utilities_kokkos<Scalar, LocalOrdinal, GlobalOrdinal, Node>::GetMatrixDiagonalInverse(
const Matrix& A, Magnitude tol) {
139 RCP<const Map> rowMap = A.getRowMap();
140 RCP<Vector> diag = VectorFactory::Build(rowMap);
141 ArrayRCP<SC> diagVals = diag->getDataNonConst(0);
143 size_t numRows = rowMap->getNodeNumElements();
145 Teuchos::ArrayView<const LO> cols;
146 Teuchos::ArrayView<const SC> vals;
147 for (
size_t i = 0; i < numRows; ++i) {
148 A.getLocalRowView(i, cols, vals);
151 for (; j < cols.size(); ++j) {
152 if (Teuchos::as<size_t>(cols[j]) == i) {
153 if(Teuchos::ScalarTraits<SC>::magnitude(vals[j]) > tol)
154 diagVals[i] = Teuchos::ScalarTraits<SC>::one() / vals[j];
156 diagVals[i]=Teuchos::ScalarTraits<SC>::zero();
160 if (j == cols.size()) {
162 diagVals[i]=Teuchos::ScalarTraits<SC>::zero();
170 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
171 Teuchos::ArrayRCP<Scalar> Utilities_kokkos<Scalar, LocalOrdinal, GlobalOrdinal, Node>::GetLumpedMatrixDiagonal(
const Matrix &A) {
173 size_t numRows = A.getRowMap()->getNodeNumElements();
174 Teuchos::ArrayRCP<SC> diag(numRows);
176 Teuchos::ArrayView<const LO> cols;
177 Teuchos::ArrayView<const SC> vals;
178 for (
size_t i = 0; i < numRows; ++i) {
179 A.getLocalRowView(i, cols, vals);
181 diag[i] = Teuchos::ScalarTraits<Scalar>::zero();
182 for (LO j = 0; j < cols.size(); ++j) {
183 diag[i] += Teuchos::ScalarTraits<Scalar>::magnitude(vals[j]);
190 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
191 RCP<Xpetra::Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Utilities_kokkos<Scalar, LocalOrdinal, GlobalOrdinal, Node>::GetMatrixOverlappedDiagonal(
const Matrix& A) {
193 RCP<const Map> rowMap = A.getRowMap(), colMap = A.getColMap();
194 RCP<Vector> localDiag = VectorFactory::Build(rowMap);
197 const CrsMatrixWrap* crsOp =
dynamic_cast<const CrsMatrixWrap*
>(&A);
199 throw Exceptions::RuntimeError(
"cast to CrsMatrixWrap failed");
201 Teuchos::ArrayRCP<size_t> offsets;
202 crsOp->getLocalDiagOffsets(offsets);
203 crsOp->getLocalDiagCopy(*localDiag,offsets());
206 ArrayRCP<SC> localDiagVals = localDiag->getDataNonConst(0);
207 Teuchos::ArrayRCP<SC> diagVals = GetMatrixDiagonal(A);
208 for (LO i = 0; i < localDiagVals.size(); i++)
209 localDiagVals[i] = diagVals[i];
210 localDiagVals = diagVals = null;
213 RCP<Vector> diagonal = VectorFactory::Build(colMap);
214 RCP< const Import> importer;
215 importer = A.getCrsGraph()->getImporter();
216 if (importer == Teuchos::null) {
217 importer = ImportFactory::Build(rowMap, colMap);
219 diagonal->doImport(*localDiag, *(importer), Xpetra::INSERT);
224 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
225 void Utilities_kokkos<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MyOldScaleMatrix(Matrix& Op,
const Teuchos::ArrayRCP<const SC>& scalingVector,
bool doInverse,
227 bool doOptimizeStorage)
229 SC one = Teuchos::ScalarTraits<SC>::one();
230 Teuchos::ArrayRCP<SC> sv(scalingVector.size());
232 for (
int i = 0; i < scalingVector.size(); ++i)
233 sv[i] = one / scalingVector[i];
235 for (
int i = 0; i < scalingVector.size(); ++i)
236 sv[i] = scalingVector[i];
239 switch (Op.getRowMap()->lib()) {
240 case Xpetra::UseTpetra:
241 MyOldScaleMatrix_Tpetra(Op, sv, doFillComplete, doOptimizeStorage);
244 case Xpetra::UseEpetra:
247 throw std::runtime_error(
"FIXME");
248 #ifndef __NVCC__ //prevent nvcc warning 253 throw Exceptions::RuntimeError(
"Only Epetra and Tpetra matrices can be scaled.");
254 #ifndef __NVCC__ //prevent nvcc warning 260 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
261 void Utilities_kokkos<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MyOldScaleMatrix_Epetra(Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Op,
const Teuchos::ArrayRCP<Scalar>& scalingVector,
bool doFillComplete,
bool doOptimizeStorage) {
262 throw Exceptions::RuntimeError(
"MyOldScaleMatrix_Epetra: Epetra needs SC=double and LO=GO=int.");
265 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
266 void Utilities_kokkos<Scalar, LocalOrdinal, GlobalOrdinal, Node>::MyOldScaleMatrix_Tpetra(Matrix& Op,
const Teuchos::ArrayRCP<SC>& scalingVector,
268 bool doOptimizeStorage)
270 #ifdef HAVE_MUELU_TPETRA 272 Tpetra::CrsMatrix<SC,LO,GO,NO>& tpOp = Op2NonConstTpetraCrs(Op);
274 const RCP<const Tpetra::Map<LO,GO,NO> > rowMap = tpOp.getRowMap();
275 const RCP<const Tpetra::Map<LO,GO,NO> > domainMap = tpOp.getDomainMap();
276 const RCP<const Tpetra::Map<LO,GO,NO> > rangeMap = tpOp.getRangeMap();
278 size_t maxRowSize = tpOp.getNodeMaxNumRowEntries();
279 if (maxRowSize == Teuchos::as<size_t>(-1))
282 std::vector<SC> scaledVals(maxRowSize);
283 if (tpOp.isFillComplete())
286 if (Op.isLocallyIndexed() ==
true) {
287 Teuchos::ArrayView<const LO> cols;
288 Teuchos::ArrayView<const SC> vals;
290 for (
size_t i = 0; i < rowMap->getNodeNumElements(); ++i) {
291 tpOp.getLocalRowView(i, cols, vals);
292 size_t nnz = tpOp.getNumEntriesInLocalRow(i);
293 if (nnz > maxRowSize) {
295 scaledVals.resize(maxRowSize);
297 for (
size_t j = 0; j < nnz; ++j)
298 scaledVals[j] = vals[j]*scalingVector[i];
301 Teuchos::ArrayView<const SC> valview(&scaledVals[0], nnz);
302 tpOp.replaceLocalValues(i, cols, valview);
307 Teuchos::ArrayView<const GO> cols;
308 Teuchos::ArrayView<const SC> vals;
310 for (
size_t i = 0; i < rowMap->getNodeNumElements(); ++i) {
311 GO gid = rowMap->getGlobalElement(i);
312 tpOp.getGlobalRowView(gid, cols, vals);
313 size_t nnz = tpOp.getNumEntriesInGlobalRow(gid);
314 if (nnz > maxRowSize) {
316 scaledVals.resize(maxRowSize);
319 for (
size_t j = 0; j < nnz; ++j)
320 scaledVals[j] = vals[j]*scalingVector[i];
323 Teuchos::ArrayView<const SC> valview(&scaledVals[0], nnz);
324 tpOp.replaceGlobalValues(gid, cols, valview);
329 if (doFillComplete) {
330 if (domainMap == Teuchos::null || rangeMap == Teuchos::null)
331 throw Exceptions::RuntimeError(
"In Utils_kokkos::Scaling: cannot fillComplete because the domain and/or range map hasn't been defined");
333 RCP<Teuchos::ParameterList> params = rcp(
new Teuchos::ParameterList());
334 params->set(
"Optimize Storage", doOptimizeStorage);
335 params->set(
"No Nonlocal Changes",
true);
336 Op.fillComplete(Op.getDomainMap(), Op.getRangeMap(), params);
339 throw Exceptions::RuntimeError(
"Only Tpetra::CrsMatrix types can be scaled (Err.1)");
342 throw Exceptions::RuntimeError(
"Matrix scaling is not possible because Tpetra has not been enabled.");
347 template <
class SC,
class LO,
class GO,
class NO>
348 Kokkos::View<const bool*, typename NO::device_type>
350 const typename Teuchos::ScalarTraits<SC>::magnitudeType& tol,
351 const bool count_twos_as_dirichlet) {
352 using ATS = Kokkos::ArithTraits<SC>;
353 using range_type = Kokkos::RangePolicy<LO, typename NO::execution_space>;
355 auto localMatrix = A.getLocalMatrix();
356 LO numRows = A.getNodeNumRows();
358 Kokkos::View<bool*, typename NO::device_type> boundaryNodes(Kokkos::ViewAllocateWithoutInitializing(
"boundaryNodes"), numRows);
359 if (count_twos_as_dirichlet)
360 Kokkos::parallel_for(
"MueLu:Utils::DetectDirichletRows_Twos_As_Dirichlet", range_type(0,numRows),
361 KOKKOS_LAMBDA(
const LO row) {
362 auto rowView = localMatrix.row(row);
363 auto length = rowView.length;
365 boundaryNodes(row) =
true;
367 decltype(length) colID;
368 for (colID = 0; colID < length; colID++)
369 if ((rowView.colidx(colID) != row) &&
370 (ATS::magnitude(rowView.value(colID)) > tol)) {
371 if (!boundaryNodes(row))
373 boundaryNodes(row) =
false;
376 boundaryNodes(row) =
true;
380 Kokkos::parallel_for(
"MueLu:Utils::DetectDirichletRows", range_type(0,numRows),
381 KOKKOS_LAMBDA(
const LO row) {
382 auto rowView = localMatrix.row(row);
383 auto length = rowView.length;
385 boundaryNodes(row) =
true;
386 for (decltype(length) colID = 0; colID < length; colID++)
387 if ((rowView.colidx(colID) != row) &&
388 (ATS::magnitude(rowView.value(colID)) > tol)) {
389 boundaryNodes(row) =
false;
394 return boundaryNodes;
397 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
398 Kokkos::View<const bool*, typename Node::device_type>
400 DetectDirichletRows(
const Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>& A,
const typename Teuchos::ScalarTraits<Scalar>::magnitudeType& tol,
const bool count_twos_as_dirichlet) {
401 return MueLu::DetectDirichletRows<Scalar,LocalOrdinal,GlobalOrdinal,Node>(A, tol, count_twos_as_dirichlet);
404 template <
class Node>
405 Kokkos::View<const bool*, typename Node::device_type>
407 DetectDirichletRows(
const Xpetra::Matrix<double,int,int,Node>& A,
const typename Teuchos::ScalarTraits<double>::magnitudeType& tol,
const bool count_twos_as_dirichlet) {
408 return MueLu::DetectDirichletRows<double,int,int,Node>(A, tol,count_twos_as_dirichlet);
412 template <
class SC,
class LO,
class GO,
class NO>
413 Kokkos::View<const bool*, typename NO::device_type>
415 const Kokkos::View<const bool*, typename NO::device_type>& dirichletRows) {
416 using ATS = Kokkos::ArithTraits<SC>;
417 using range_type = Kokkos::RangePolicy<LO, typename NO::execution_space>;
419 SC zero = ATS::zero();
422 auto localMatrix = A.getLocalMatrix();
423 LO numRows = A.getNodeNumRows();
425 Teuchos::RCP<const Xpetra::Map<LO,GO,NO> > domMap = A.getDomainMap();
426 Teuchos::RCP<const Xpetra::Map<LO,GO,NO> > colMap = A.getColMap();
427 Teuchos::RCP<Xpetra::MultiVector<SC,LO,GO,NO> > myColsToZero = Xpetra::MultiVectorFactory<SC,LO,GO,NO>::Build(colMap,1);
428 myColsToZero->putScalar(zero);
429 auto myColsToZeroView = myColsToZero->template getLocalView<typename NO::device_type>();
431 Kokkos::parallel_for(
"MueLu:Utils::DetectDirichletCols1", range_type(0,numRows),
432 KOKKOS_LAMBDA(
const LO row) {
433 if (dirichletRows(row)) {
434 auto rowView = localMatrix.row(row);
435 auto length = rowView.length;
437 for (decltype(length) colID = 0; colID < length; colID++)
438 myColsToZeroView(rowView.colidx(colID),0) = one;
442 Teuchos::RCP<Xpetra::MultiVector<SC,LO,GO,NO> > globalColsToZero = Xpetra::MultiVectorFactory<SC,LO,GO,NO>::Build(domMap,1);
443 globalColsToZero->putScalar(zero);
444 Teuchos::RCP<Xpetra::Export<LO,GO,NO> > exporter = Xpetra::ExportFactory<LO,GO,NO>::Build(colMap,domMap);
446 globalColsToZero->doExport(*myColsToZero,*exporter,Xpetra::ADD);
448 myColsToZero->doImport(*globalColsToZero,*exporter,Xpetra::INSERT);
450 auto myCols = myColsToZero->template getLocalView<typename NO::device_type>();
451 size_t numColEntries = colMap->getNodeNumElements();
452 Kokkos::View<bool*, typename NO::device_type> dirichletCols(Kokkos::ViewAllocateWithoutInitializing(
"dirichletCols"), numColEntries);
453 const typename ATS::magnitudeType eps = 2.0*ATS::eps();
455 Kokkos::parallel_for(
"MueLu:Utils::DetectDirichletCols2", range_type(0,numColEntries),
456 KOKKOS_LAMBDA(
const size_t i) {
457 dirichletCols(i) = ATS::magnitude(myCols(i,0))>eps;
459 return dirichletCols;
463 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
464 Kokkos::View<const bool*, typename Node::device_type>
467 const Kokkos::View<const bool*, typename NO::device_type>& dirichletRows) {
468 return MueLu::DetectDirichletCols<Scalar,LocalOrdinal,GlobalOrdinal,Node>(A, dirichletRows);
471 template <
class Node>
472 Kokkos::View<const bool*, typename Node::device_type>
475 const Kokkos::View<const bool*, typename Node::device_type>& dirichletRows) {
476 return MueLu::DetectDirichletCols<double,int,int,Node>(A, dirichletRows);
481 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
484 const Kokkos::View<const bool*, typename Node::device_type>& dirichletRows,
485 Scalar replaceWith) {
486 using ATS = Kokkos::ArithTraits<Scalar>;
487 using range_type = Kokkos::RangePolicy<LocalOrdinal, typename Node::execution_space>;
489 auto localMatrix = A->getLocalMatrix();
490 LocalOrdinal numRows = A->getNodeNumRows();
492 Kokkos::parallel_for(
"MueLu:Utils::ZeroDirichletRows", range_type(0,numRows),
493 KOKKOS_LAMBDA(
const LocalOrdinal row) {
494 if (dirichletRows(row)) {
495 auto rowView = localMatrix.row(row);
496 auto length = rowView.length;
497 for (decltype(length) colID = 0; colID < length; colID++)
498 rowView.value(colID) = replaceWith;
503 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
506 ZeroDirichletRows(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& A,
507 const Kokkos::View<const bool*, typename NO::device_type>& dirichletRows,
508 Scalar replaceWith) {
509 MueLu::ZeroDirichletRows<Scalar,LocalOrdinal,GlobalOrdinal,Node>(A, dirichletRows, replaceWith);
512 template <
class Node>
516 const Kokkos::View<const bool*, typename Node::device_type>& dirichletRows,
517 double replaceWith) {
518 return MueLu::ZeroDirichletRows<double,int,int,Node>(A, dirichletRows, replaceWith);
523 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
526 const Kokkos::View<const bool*, typename Node::device_type>& dirichletRows,
527 Scalar replaceWith) {
528 using range_type = Kokkos::RangePolicy<LocalOrdinal, typename Node::execution_space>;
529 auto myCols = X->template getLocalView<typename Node::device_type>();
530 size_t numVecs = X->getNumVectors();
531 Kokkos::parallel_for(
"MueLu:Utils::ZeroDirichletRows_MV", range_type(0,dirichletRows.size()),
532 KOKKOS_LAMBDA(
const size_t i) {
533 if (dirichletRows(i)) {
534 for(
size_t j=0; j<numVecs; j++)
535 myCols(i,j) = replaceWith;
540 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
543 ZeroDirichletRows(RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> >& X,
544 const Kokkos::View<const bool*, typename NO::device_type>& dirichletRows,
545 Scalar replaceWith) {
546 MueLu::ZeroDirichletRows<Scalar,LocalOrdinal,GlobalOrdinal,Node>(X, dirichletRows, replaceWith);
549 template <
class Node>
553 const Kokkos::View<const bool*, typename Node::device_type>& dirichletRows,
554 double replaceWith) {
555 return MueLu::ZeroDirichletRows<double,int,int,Node>(X, dirichletRows, replaceWith);
560 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
563 const Kokkos::View<const bool*, typename Node::device_type>& dirichletCols,
564 Scalar replaceWith) {
565 using ATS = Kokkos::ArithTraits<Scalar>;
566 using range_type = Kokkos::RangePolicy<LocalOrdinal, typename Node::execution_space>;
568 auto localMatrix = A->getLocalMatrix();
569 LocalOrdinal numRows = A->getNodeNumRows();
571 Kokkos::parallel_for(
"MueLu:Utils::ZeroDirichletCols", range_type(0,numRows),
572 KOKKOS_LAMBDA(
const LocalOrdinal row) {
573 auto rowView = localMatrix.row(row);
574 auto length = rowView.length;
575 for (decltype(length) colID = 0; colID < length; colID++)
576 if (dirichletCols(rowView.colidx(colID))) {
577 rowView.value(colID) = replaceWith;
582 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
585 ZeroDirichletCols(RCP<Xpetra::Matrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& A,
586 const Kokkos::View<const bool*, typename NO::device_type>& dirichletCols,
587 Scalar replaceWith) {
588 MueLu::ZeroDirichletCols<Scalar,LocalOrdinal,GlobalOrdinal,Node>(A, dirichletCols, replaceWith);
591 template <
class Node>
595 const Kokkos::View<const bool*, typename Node::device_type>& dirichletCols,
596 double replaceWith) {
597 return MueLu::ZeroDirichletCols<double,int,int,Node>(A, dirichletCols, replaceWith);
601 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
602 RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
603 Utilities_kokkos<Scalar, LocalOrdinal, GlobalOrdinal, Node>::
604 RealValuedToScalarMultiVector(RCP<RealValuedMultiVector > X) {
605 RCP<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > Xscalar;
606 #if defined(HAVE_XPETRA_TPETRA) && defined(HAVE_TPETRA_INST_COMPLEX_DOUBLE) 607 using range_type = Kokkos::RangePolicy<LocalOrdinal, typename Node::execution_space>;
610 if (
typeid(Scalar).name() ==
typeid(std::complex<double>).name()) {
611 size_t numVecs = X->getNumVectors();
612 Xscalar = Xpetra::MultiVectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node>::Build(X->getMap(),numVecs);
613 auto XVec = X->template getLocalView<typename Node::device_type>();
614 auto XVecScalar = Xscalar->template getLocalView<typename Node::device_type>();
616 Kokkos::parallel_for(
"MueLu:Utils::RealValuedToScalarMultiVector", range_type(0,X->getLocalLength()),
617 KOKKOS_LAMBDA(
const size_t i) {
618 for (
size_t j=0; j<numVecs; j++)
619 XVecScalar(i,j) = XVec(i,j);
623 Xscalar = rcp_dynamic_cast<Xpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> >(X);
627 template <
class Node>
628 RCP<Xpetra::MultiVector<double,int,int,Node> >
629 Utilities_kokkos<double,int,int,Node>::
630 RealValuedToScalarMultiVector(RCP<Xpetra::MultiVector<Magnitude,int,int,Node> > X) {
636 #define MUELU_UTILITIES_KOKKOS_SHORT 637 #endif // MUELU_UTILITIES_KOKKOS_DEF_HPP
void ZeroDirichletCols(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Kokkos::View< const bool *, typename Node::device_type > &dirichletCols, Scalar replaceWith)
void ZeroDirichletRows(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Kokkos::View< const bool *, typename Node::device_type > &dirichletRows, Scalar replaceWith)
Namespace for MueLu classes and methods.
void ZeroDirichletRows(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &X, const Kokkos::View< const bool *, typename Node::device_type > &dirichletRows, Scalar replaceWith)
Kokkos::View< const bool *, typename NO::device_type > DetectDirichletCols(const Xpetra::Matrix< SC, LO, GO, NO > &A, const Kokkos::View< const bool *, typename NO::device_type > &dirichletRows)
Kokkos::View< const bool *, typename NO::device_type > DetectDirichletRows(const Xpetra::Matrix< SC, LO, GO, NO > &A, const typename Teuchos::ScalarTraits< SC >::magnitudeType &tol, const bool count_twos_as_dirichlet)