46 #ifndef MUELU_REFMAXWELL_DEF_HPP 47 #define MUELU_REFMAXWELL_DEF_HPP 59 #include "MueLu_AmalgamationFactory.hpp" 60 #include "MueLu_RAPFactory.hpp" 61 #include "MueLu_ThresholdAFilterFactory.hpp" 62 #include "MueLu_TransPFactory.hpp" 63 #include "MueLu_SmootherFactory.hpp" 65 #include "MueLu_CoalesceDropFactory.hpp" 66 #include "MueLu_CoarseMapFactory.hpp" 67 #include "MueLu_CoordinatesTransferFactory.hpp" 68 #include "MueLu_UncoupledAggregationFactory.hpp" 69 #include "MueLu_TentativePFactory.hpp" 70 #include "MueLu_Utilities.hpp" 72 #ifdef HAVE_MUELU_KOKKOS_REFACTOR 73 #include "MueLu_CoalesceDropFactory_kokkos.hpp" 74 #include "MueLu_CoarseMapFactory_kokkos.hpp" 75 #include "MueLu_CoordinatesTransferFactory_kokkos.hpp" 76 #include "MueLu_UncoupledAggregationFactory_kokkos.hpp" 77 #include "MueLu_TentativePFactory_kokkos.hpp" 78 #include "MueLu_Utilities_kokkos.hpp" 79 #include <Kokkos_Core.hpp> 80 #include <KokkosSparse_CrsMatrix.hpp> 85 #include "MueLu_MLParameterListInterpreter.hpp" 86 #include "MueLu_ParameterListInterpreter.hpp" 91 #ifdef HAVE_MUELU_CUDA 92 #include "cuda_profiler_api.h" 98 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
100 return SM_Matrix_->getDomainMap();
104 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
106 return SM_Matrix_->getRangeMap();
110 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
113 parameterList_ = list;
114 disable_addon_ = list.
get(
"refmaxwell: disable addon",
true);
115 mode_ = list.
get(
"refmaxwell: mode",
"additive");
116 use_as_preconditioner_ = list.
get<
bool>(
"refmaxwell: use as preconditioner");
117 dump_matrices_ = list.
get(
"refmaxwell: dump matrices",
false);
120 precList11_ = list.
sublist(
"refmaxwell: 11list");
123 precList22_ = list.
sublist(
"refmaxwell: 22list");
126 smootherList_ = list.
sublist(
"smoother: params");
129 #if !defined(HAVE_MUELU_KOKKOS_REFACTOR) 131 #elif defined(HAVE_MUELU_KOKKOS_REFACTOR_USE_BY_DEFAULT) 132 useKokkos_ = list.
get(
"use kokkos refactor",
true);
134 useKokkos_ = list.
get(
"use kokkos refactor",
false);
140 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
144 #ifdef HAVE_MUELU_CUDA 145 if (parameterList_.get<
bool>(
"refmaxwell: cuda profile setup",
false)) cudaProfilerStart();
150 bool defaultFilter =
false;
155 if (parameterList_.get<
bool>(
"refmaxwell: filter D0",
true) && D0_Matrix_->getNodeMaxNumRowEntries()>2) {
159 fineLevel.
Set(
"A",D0_Matrix_);
160 fineLevel.
setlib(D0_Matrix_->getDomainMap()->lib());
164 ThreshFact->
Build(fineLevel);
168 defaultFilter =
true;
171 if (parameterList_.get<
bool>(
"refmaxwell: filter SM", defaultFilter)) {
172 RCP<Vector> diag = VectorFactory::Build(SM_Matrix_->getRowMap());
174 SM_Matrix_->getLocalDiagCopy(*diag);
180 fineLevel.
Set(
"A",SM_Matrix_);
181 fineLevel.
setlib(SM_Matrix_->getDomainMap()->lib());
184 ThreshFact->
Build(fineLevel);
188 if (parameterList_.get<
bool>(
"refmaxwell: filter M1", defaultFilter)) {
189 RCP<Vector> diag = VectorFactory::Build(M1_Matrix_->getRowMap());
191 M1_Matrix_->getLocalDiagCopy(*diag);
197 fineLevel.
Set(
"A",M1_Matrix_);
198 fineLevel.
setlib(M1_Matrix_->getDomainMap()->lib());
201 ThreshFact->
Build(fineLevel);
209 #ifdef HAVE_MUELU_KOKKOS_REFACTOR 223 if(Nullspace_ != Teuchos::null) {
226 else if(Nullspace_ == Teuchos::null && Coords_ != Teuchos::null) {
231 Coords_->norm2(norms);
232 for (
size_t i=0;i<Coords_->getNumVectors();i++)
233 norms[i] = ((realMagnitudeType)1.0)/norms[i];
234 Coords_->scale(norms());
235 Nullspace_ = MultiVectorFactory::Build(SM_Matrix_->getRowMap(),Coords_->getNumVectors());
238 #ifdef HAVE_MUELU_KOKKOS_REFACTOR 241 CoordsSC = Utilities_kokkos::RealValuedToScalarMultiVector(Coords_);
247 D0_Matrix_->apply(*CoordsSC,*Nullspace_);
250 GetOStream(
Errors) <<
"MueLu::RefMaxwell::compute(): either the nullspace or the nodal coordinates must be provided." << std::endl;
254 #ifdef HAVE_MUELU_KOKKOS_REFACTOR 267 if(P11_==Teuchos::null) {
269 Level fineLevel, coarseLevel;
275 fineLevel.
Set(
"A",M1_Matrix_);
276 coarseLevel.
Set(
"P",D0_Matrix_);
277 coarseLevel.
setlib(M1_Matrix_->getDomainMap()->lib());
278 fineLevel.
setlib(M1_Matrix_->getDomainMap()->lib());
284 rapList.
set(
"transpose: use implicit", parameterList_.get<
bool>(
"transpose: use implicit",
false));
285 rapList.
set(
"rap: fix zero diagonals", parameterList_.get<
bool>(
"rap: fix zero diagonals",
true));
286 rapList.
set(
"rap: triple product", parameterList_.get<
bool>(
"rap: triple product",
false));
290 if (!parameterList_.get<
bool>(
"transpose: use implicit",
false)) {
300 GetOStream(
Runtime0) <<
"RefMaxwell::compute(): building special prolongator" << std::endl;
303 #ifdef HAVE_MUELU_KOKKOS_REFACTOR 305 R11_ = Utilities_kokkos::Transpose(*P11_);
317 #ifdef HAVE_MUELU_KOKKOS_REFACTOR 321 if (useKokkos_ && precList11_.isParameter(
"aggregation: drop tol") && precList11_.get<
double>(
"aggregation: drop tol") != 0.0) {
322 GetOStream(
Warnings0) <<
"RefMaxwell::compute(): Setting \"aggregation: drop tol\". to 0.0, since CoalesceDropFactory_kokkos does not " 323 <<
"support BlockSize > 1 and drop tol != 0.0" << std::endl;
324 precList11_.set(
"aggregation: drop tol", 0.0);
331 GetOStream(
Runtime0) <<
"RefMaxwell::compute(): nuking BC edges of D0" << std::endl;
333 D0_Matrix_->resumeFill();
336 #ifdef HAVE_MUELU_KOKKOS_REFACTOR 348 D0_Matrix_->fillComplete(D0_Matrix_->getDomainMap(),D0_Matrix_->getRangeMap());
352 GetOStream(
Runtime0) <<
"RefMaxwell::compute(): building MG for (2,2)-block" << std::endl;
357 Level fineLevel, coarseLevel;
363 fineLevel.
Set(
"A",SM_Matrix_);
364 coarseLevel.
Set(
"P",D0_Matrix_);
365 coarseLevel.
setlib(SM_Matrix_->getDomainMap()->lib());
366 fineLevel.
setlib(SM_Matrix_->getDomainMap()->lib());
372 rapList.
set(
"transpose: use implicit",
false);
373 rapList.
set(
"rap: fix zero diagonals", parameterList_.get<
bool>(
"rap: fix zero diagonals",
true));
374 rapList.
set(
"rap: triple product", parameterList_.get<
bool>(
"rap: triple product",
false));
381 coarseLevel.
Request(
"R", transPFactory.
get());
384 A22_->setObjectLabel(
"RefMaxwell (2,2)");
392 if (parameterList_.isType<std::string>(
"smoother: type") &&
393 parameterList_.get<std::string>(
"smoother: type") ==
"hiptmair" &&
397 #if defined(HAVE_MUELU_IFPACK2) && (!defined(HAVE_MUELU_EPETRA) || (defined(HAVE_MUELU_INST_DOUBLE_INT_INT))) 400 if (smootherList_.isSublist(
"smoother: pre params"))
401 smootherPreList = smootherList_.
sublist(
"smoother: pre params");
402 else if (smootherList_.isSublist(
"smoother: params"))
403 smootherPreList = smootherList_.
sublist(
"smoother: params");
404 hiptmairPreList.
set(
"hiptmair: smoother type 1",
405 smootherPreList.
get<std::string>(
"hiptmair: smoother type 1",
"CHEBYSHEV"));
406 hiptmairPreList.
set(
"hiptmair: smoother type 2",
407 smootherPreList.
get<std::string>(
"hiptmair: smoother type 2",
"CHEBYSHEV"));
408 if(smootherPreList.
isSublist(
"hiptmair: smoother list 1"))
409 hiptmairPreList.
set(
"hiptmair: smoother list 1", smootherPreList.
sublist(
"hiptmair: smoother list 1"));
410 if(smootherPreList.
isSublist(
"hiptmair: smoother list 2"))
411 hiptmairPreList.
set(
"hiptmair: smoother list 2", smootherPreList.
sublist(
"hiptmair: smoother list 2"));
412 hiptmairPreList.
set(
"hiptmair: pre or post",
413 smootherPreList.
get<std::string>(
"hiptmair: pre or post",
"pre"));
414 hiptmairPreList.
set(
"hiptmair: zero starting solution",
415 smootherPreList.
get<
bool>(
"hiptmair: zero starting solution",
true));
417 if (smootherList_.isSublist(
"smoother: post params"))
418 smootherPostList = smootherList_.
sublist(
"smoother: post params");
419 else if (smootherList_.isSublist(
"smoother: params"))
420 smootherPostList = smootherList_.
sublist(
"smoother: params");
421 hiptmairPostList.
set(
"hiptmair: smoother type 1",
422 smootherPostList.
get<std::string>(
"hiptmair: smoother type 1",
"CHEBYSHEV"));
423 hiptmairPostList.
set(
"hiptmair: smoother type 2",
424 smootherPostList.
get<std::string>(
"hiptmair: smoother type 2",
"CHEBYSHEV"));
425 if(smootherPostList.
isSublist(
"hiptmair: smoother list 1"))
426 hiptmairPostList.
set(
"hiptmair: smoother list 1", smootherPostList.
sublist(
"hiptmair: smoother list 1"));
427 if(smootherPostList.
isSublist(
"hiptmair: smoother list 2"))
428 hiptmairPostList.
set(
"hiptmair: smoother list 2", smootherPostList.
sublist(
"hiptmair: smoother list 2"));
429 hiptmairPostList.
set(
"hiptmair: pre or post",
430 smootherPostList.
get<std::string>(
"hiptmair: pre or post",
"post"));
431 hiptmairPostList.
set(
"hiptmair: zero starting solution",
432 smootherPostList.
get<
bool>(
"hiptmair: zero starting solution",
false));
434 typedef Tpetra::RowMatrix<SC, LO, GO, NO> TROW;
439 hiptmairPreSmoother_ =
Teuchos::rcp(
new Ifpack2::Hiptmair<TROW>(EdgeMatrix,NodeMatrix,PMatrix) );
440 hiptmairPreSmoother_ -> setParameters(hiptmairPreList);
441 hiptmairPreSmoother_ -> initialize();
442 hiptmairPreSmoother_ -> compute();
443 hiptmairPostSmoother_ =
Teuchos::rcp(
new Ifpack2::Hiptmair<TROW>(EdgeMatrix,NodeMatrix,PMatrix) );
444 hiptmairPostSmoother_ -> setParameters(hiptmairPostList);
445 hiptmairPostSmoother_ -> initialize();
446 hiptmairPostSmoother_ -> compute();
447 useHiptmairSmoothing_ =
true;
450 #endif // defined(HAVE_MUELU_IFPACK2) && (!defined(HAVE_MUELU_EPETRA) || defined(HAVE_MUELU_INST_DOUBLE_INT_INT)) 452 if (parameterList_.isType<std::string>(
"smoother: pre type") && parameterList_.isType<std::string>(
"smoother: post type")) {
453 std::string preSmootherType = parameterList_.get<std::string>(
"smoother: pre type");
454 std::string postSmootherType = parameterList_.get<std::string>(
"smoother: post type");
457 if (parameterList_.isSublist(
"smoother: pre params"))
458 preSmootherList = parameterList_.
sublist(
"smoother: pre params");
459 if (parameterList_.isSublist(
"smoother: post params"))
460 postSmootherList = parameterList_.
sublist(
"smoother: post params");
467 level.
Set(
"A",SM_Matrix_);
468 level.
setlib(SM_Matrix_->getDomainMap()->lib());
476 level.
Request(
"PreSmoother",preSmootherFact.
get());
477 preSmootherFact->
Build(level);
480 level.
Request(
"PostSmoother",postSmootherFact.
get());
481 postSmootherFact->
Build(level);
484 std::string smootherType = parameterList_.
get<std::string>(
"smoother: type",
"CHEBYSHEV");
487 level.SetFactoryManager(factoryHandler);
489 level.setObjectLabel(
"RefMaxwell (1,1)");
490 level.Set(
"A",SM_Matrix_);
491 level.setlib(SM_Matrix_->getDomainMap()->lib());
494 level.Request(
"PreSmoother",SmootherFact.
get());
495 SmootherFact->
Build(level);
497 PostSmoother_ = PreSmoother_;
499 useHiptmairSmoothing_ =
false;
504 P11res_ = MultiVectorFactory::Build(P11_->getDomainMap(),1);
505 P11x_ = MultiVectorFactory::Build(P11_->getDomainMap(),1);
506 D0res_ = MultiVectorFactory::Build(D0_Matrix_->getDomainMap(),1);
507 D0x_ = MultiVectorFactory::Build(D0_Matrix_->getDomainMap(),1);
508 residual_ = MultiVectorFactory::Build(SM_Matrix_->getDomainMap(),1);
510 #ifdef HAVE_MUELU_CUDA 511 if (parameterList_.get<
bool>(
"refmaxwell: cuda profile setup",
false)) cudaProfilerStop();
516 if (dump_matrices_) {
517 GetOStream(
Runtime0) <<
"RefMaxwell::compute(): dumping data" << std::endl;
521 #ifndef HAVE_MUELU_KOKKOS_REFACTOR 522 std::ofstream outBCrows(
"BCrows.mat");
523 std::copy(BCrows_.begin(), BCrows_.end(), std::ostream_iterator<LO>(outBCrows,
"\n"));
524 std::ofstream outBCcols(
"BCcols.mat");
525 std::copy(BCcols_.begin(), BCcols_.end(), std::ostream_iterator<LO>(outBCcols,
"\n"));
528 std::ofstream outBCrows(
"BCrows.mat");
529 auto BCrows = Kokkos::create_mirror_view (BCrowsKokkos_);
531 for (
size_t i = 0; i < BCrows.size(); i++)
532 outBCrows << BCrows[i] <<
"\n";
534 std::ofstream outBCcols(
"BCcols.mat");
535 auto BCcols = Kokkos::create_mirror_view (BCcolsKokkos_);
537 for (
size_t i = 0; i < BCcols.size(); i++)
538 outBCcols << BCcols[i] <<
"\n";
540 std::ofstream outBCrows(
"BCrows.mat");
541 std::copy(BCrows_.begin(), BCrows_.end(), std::ostream_iterator<LO>(outBCrows,
"\n"));
542 std::ofstream outBCcols(
"BCcols.mat");
543 std::copy(BCcols_.begin(), BCcols_.end(), std::ostream_iterator<LO>(outBCcols,
"\n"));
547 if (Coords_ != Teuchos::null)
558 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
573 size_t dim = Nullspace_->getNumVectors();
574 size_t numLocalRows = SM_Matrix_->getNodeNumRows();
579 bool read_P_from_file = parameterList_.
get(
"refmaxwell: read_P_from_file",
false);
580 if (read_P_from_file) {
583 std::string P_filename = parameterList_.get(
"refmaxwell: P_filename",std::string(
"P.mat"));
586 Level fineLevel, coarseLevel;
592 fineLevel.
Set(
"A",A_nodal_Matrix_);
593 fineLevel.
Set(
"Coordinates",Coords_);
594 fineLevel.
Set(
"DofsPerNode",1);
595 coarseLevel.
setlib(A_nodal_Matrix_->getDomainMap()->lib());
596 fineLevel.
setlib(A_nodal_Matrix_->getDomainMap()->lib());
600 LocalOrdinal NSdim = 1;
601 RCP<MultiVector> nullSpace = MultiVectorFactory::Build(A_nodal_Matrix_->getRowMap(),NSdim);
602 nullSpace->putScalar(SC_ONE);
603 fineLevel.
Set(
"Nullspace",nullSpace);
606 #ifdef HAVE_MUELU_KOKKOS_REFACTOR 607 RCP<Factory> dropFact, UncoupledAggFact, coarseMapFact, TentativePFact, Tfact;
609 dropFact =
rcp(
new CoalesceDropFactory_kokkos());
610 UncoupledAggFact =
rcp(
new UncoupledAggregationFactory_kokkos());
611 coarseMapFact =
rcp(
new CoarseMapFactory_kokkos());
612 TentativePFact =
rcp(
new TentativePFactory_kokkos());
613 Tfact =
rcp(
new CoordinatesTransferFactory_kokkos());
628 dropFact->
SetFactory(
"UnAmalgamationInfo", amalgFact);
629 double dropTol = parameterList_.get(
"aggregation: drop tol",0.0);
632 UncoupledAggFact->
SetFactory(
"Graph", dropFact);
634 coarseMapFact->
SetFactory(
"Aggregates", UncoupledAggFact);
636 TentativePFact->
SetFactory(
"Aggregates", UncoupledAggFact);
637 TentativePFact->
SetFactory(
"UnAmalgamationInfo", amalgFact);
638 TentativePFact->
SetFactory(
"CoarseMap", coarseMapFact);
640 Tfact->
SetFactory(
"Aggregates", UncoupledAggFact);
641 Tfact->
SetFactory(
"CoarseMap", coarseMapFact);
643 coarseLevel.
Request(
"P",TentativePFact.
get());
644 coarseLevel.
Request(
"Coordinates",Tfact.
get());
646 coarseLevel.
Get(
"P",P_nodal,TentativePFact.
get());
647 coarseLevel.
Get(
"Coordinates",CoordsH_,Tfact.
get());
652 RCP<CrsMatrix> D0Crs = rcp_dynamic_cast<CrsMatrixWrap>(D0_Matrix_)->getCrsMatrix();
656 int numProcs = P_nodal->getDomainMap()->getComm()->getSize();
660 P_nodal_temp =
Teuchos::rcp(
new CrsMatrixWrap(targetMap, 0));
663 P_nodal_temp->fillComplete(rcp_dynamic_cast<CrsMatrixWrap>(P_nodal)->getCrsMatrix()->getDomainMap(),
664 rcp_dynamic_cast<CrsMatrixWrap>(P_nodal)->getCrsMatrix()->getRangeMap());
665 P_nodal_imported = P_nodal_temp->getCrsMatrix();
669 P_nodal_imported = rcp_dynamic_cast<CrsMatrixWrap>(P_nodal)->getCrsMatrix();
671 #ifdef HAVE_MUELU_KOKKOS_REFACTOR 674 using ATS = Kokkos::ArithTraits<SC>;
677 typedef typename Matrix::local_matrix_type KCRS;
678 typedef typename KCRS::device_type device_t;
679 typedef typename KCRS::StaticCrsGraphType graph_t;
680 typedef typename graph_t::row_map_type::non_const_type lno_view_t;
681 typedef typename graph_t::entries_type::non_const_type lno_nnz_view_t;
682 typedef typename KCRS::values_type::non_const_type scalar_view_t;
685 auto localP = P_nodal_imported->getLocalMatrix();
686 auto localD0 = D0_Matrix_->getLocalMatrix();
691 std::string defaultAlgo =
"mat-mat";
692 std::string algo = parameterList_.
get(
"refmaxwell: prolongator compute algorithm",defaultAlgo);
694 if (algo ==
"mat-mat") {
699 auto localD0P = D0_P_nodal->getLocalMatrix();
705 lno_view_t P11rowptr(
"P11_rowptr", numLocalRows+1);
706 lno_nnz_view_t P11colind(
"P11_colind",dim*localD0P.graph.entries.size());
707 scalar_view_t P11vals(
"P11_vals",dim*localD0P.graph.entries.size());
710 Kokkos::parallel_for(
"MueLu:RefMaxwell::buildProlongator_adjustRowptr", range_type(0,numLocalRows+1),
711 KOKKOS_LAMBDA(
const size_t i) {
712 P11rowptr(i) = dim*localD0P.graph.row_map(i);
716 Kokkos::parallel_for(
"MueLu:RefMaxwell::buildProlongator_adjustColind", range_type(0,localD0P.graph.entries.size()),
717 KOKKOS_LAMBDA(
const size_t jj) {
718 for (
size_t k = 0; k < dim; k++) {
719 P11colind(dim*jj+k) = dim*localD0P.graph.entries(jj)+k;
720 P11vals(dim*jj+k) = SC_ZERO;
724 auto localNullspace = Nullspace_->template getLocalView<device_t>();
727 if (D0_Matrix_->getNodeMaxNumRowEntries()>2) {
731 GetOStream(
Warnings0) <<
"RefMaxwell::buildProlongator(): D0 matrix has more than 2 entries per row. Taking inefficient code path." << std::endl;
735 Kokkos::parallel_for(
"MueLu:RefMaxwell::buildProlongator_enterValues_D0wZeros", range_type(0,numLocalRows),
736 KOKKOS_LAMBDA(
const size_t i) {
737 for (
size_t ll = localD0.graph.row_map(i); ll < localD0.graph.row_map(i+1); ll++) {
738 LO l = localD0.graph.entries(ll);
739 SC p = localD0.values(ll);
740 if (ATS::magnitude(p) < tol)
742 for (
size_t jj = localP.graph.row_map(l); jj < localP.graph.row_map(l+1); jj++) {
743 LO j = localP.graph.entries(jj);
744 SC v = localP.values(jj);
745 for (
size_t k = 0; k < dim; k++) {
747 SC n = localNullspace(i,k);
749 for (m = P11rowptr(i); m < P11rowptr(i+1); m++)
750 if (P11colind(m) == jNew)
752 #if defined(HAVE_MUELU_DEBUG) && !defined(HAVE_MUELU_CUDA) 755 P11vals(m) += 0.5 * v * n;
762 Kokkos::parallel_for(
"MueLu:RefMaxwell::buildProlongator_enterValues", range_type(0,numLocalRows),
763 KOKKOS_LAMBDA(
const size_t i) {
764 for (
size_t ll = localD0.graph.row_map(i); ll < localD0.graph.row_map(i+1); ll++) {
765 LO l = localD0.graph.entries(ll);
766 for (
size_t jj = localP.graph.row_map(l); jj < localP.graph.row_map(l+1); jj++) {
767 LO j = localP.graph.entries(jj);
768 SC v = localP.values(jj);
769 for (
size_t k = 0; k < dim; k++) {
771 SC n = localNullspace(i,k);
773 for (m = P11rowptr(i); m < P11rowptr(i+1); m++)
774 if (P11colind(m) == jNew)
776 #if defined(HAVE_MUELU_DEBUG) && !defined(HAVE_MUELU_CUDA) 779 P11vals(m) += 0.5 * v * n;
787 RCP<CrsMatrix> P11Crs = rcp_dynamic_cast<CrsMatrixWrap>(P11_)->getCrsMatrix();
788 P11Crs->setAllValues(P11rowptr, P11colind, P11vals);
789 P11Crs->expertStaticFillComplete(blockDomainMap, SM_Matrix_->getRangeMap());
792 TEUCHOS_TEST_FOR_EXCEPTION(
false,std::invalid_argument,algo <<
" is not a valid option for \"refmaxwell: prolongator compute algorithm\"");
794 #endif // ifdef(HAVE_MUELU_KOKKOS_REFACTOR) 799 for(
size_t i=0; i<dim; i++) {
800 nullspaceRCP[i] = Nullspace_->getData(i);
801 nullspace[i] = nullspaceRCP[i]();
812 P_nodal_imported->getAllValues(Prowptr_RCP, Pcolind_RCP, Pvals_RCP);
813 rcp_dynamic_cast<CrsMatrixWrap>(D0_Matrix_)->getCrsMatrix()->getAllValues(D0rowptr_RCP, D0colind_RCP, D0vals_RCP);
821 Prowptr = Prowptr_RCP(); Pcolind = Pcolind_RCP(); Pvals = Pvals_RCP();
822 D0rowptr = D0rowptr_RCP(); D0colind = D0colind_RCP(); D0vals = D0vals_RCP();
830 std::string defaultAlgo =
"gustavson";
831 std::string algo = parameterList_.get(
"refmaxwell: prolongator compute algorithm",defaultAlgo);
833 if (algo ==
"mat-mat") {
841 rcp_dynamic_cast<CrsMatrixWrap>(D0_P_nodal)->getCrsMatrix()->getAllValues(D0Prowptr_RCP, D0Pcolind_RCP, D0Pvals_RCP);
848 D0Prowptr = D0Prowptr_RCP(); D0Pcolind = D0Pcolind_RCP();
854 RCP<CrsMatrix> P11Crs = rcp_dynamic_cast<CrsMatrixWrap>(P11_)->getCrsMatrix();
855 P11Crs->allocateAllValues(dim*D0Pcolind.
size(), P11rowptr_RCP, P11colind_RCP, P11vals_RCP);
862 for (
size_t i = 0; i < numLocalRows+1; i++) {
863 P11rowptr[i] = dim*D0Prowptr[i];
868 for (
size_t jj = 0; jj < (size_t) D0Pcolind.
size(); jj++)
869 for (
size_t k = 0; k < dim; k++) {
870 P11colind[nnz] = dim*D0Pcolind[jj]+k;
871 P11vals[nnz] = SC_ZERO;
876 if (D0_Matrix_->getNodeMaxNumRowEntries()>2) {
880 GetOStream(
Warnings0) <<
"RefMaxwell::buildProlongator(): D0 matrix has more than 2 entries per row. Taking inefficient code path." << std::endl;
883 for (
size_t i = 0; i < numLocalRows; i++) {
884 for (
size_t ll = D0rowptr[i]; ll < D0rowptr[i+1]; ll++) {
889 for (
size_t jj = Prowptr[l]; jj < Prowptr[l+1]; jj++) {
892 for (
size_t k = 0; k < dim; k++) {
894 SC n = nullspace[k][i];
896 for (m = P11rowptr[i]; m < P11rowptr[i+1]; m++)
897 if (P11colind[m] == jNew)
899 #ifdef HAVE_MUELU_DEBUG 902 P11vals[m] += 0.5 * v * n;
909 for (
size_t i = 0; i < numLocalRows; i++) {
910 for (
size_t ll = D0rowptr[i]; ll < D0rowptr[i+1]; ll++) {
912 for (
size_t jj = Prowptr[l]; jj < Prowptr[l+1]; jj++) {
915 for (
size_t k = 0; k < dim; k++) {
917 SC n = nullspace[k][i];
919 for (m = P11rowptr[i]; m < P11rowptr[i+1]; m++)
920 if (P11colind[m] == jNew)
922 #ifdef HAVE_MUELU_DEBUG 925 P11vals[m] += 0.5 * v * n;
932 P11Crs->setAllValues(P11rowptr_RCP, P11colind_RCP, P11vals_RCP);
933 P11Crs->expertStaticFillComplete(blockDomainMap, SM_Matrix_->getRangeMap());
935 }
else if (algo ==
"gustavson") {
937 LO maxP11col = dim * P_nodal_imported->getColMap()->getMaxLocalIndex();
941 size_t nnz_alloc = dim*D0vals_RCP.
size();
947 RCP<CrsMatrix> P11Crs = rcp_dynamic_cast<CrsMatrixWrap>(P11_)->getCrsMatrix();
948 P11Crs->allocateAllValues(nnz_alloc, P11rowptr_RCP, P11colind_RCP, P11vals_RCP);
955 if (D0_Matrix_->getNodeMaxNumRowEntries()>2) {
959 GetOStream(
Warnings0) <<
"RefMaxwell::buildProlongator(): D0 matrix has more than 2 entries per row. Taking inefficient code path." << std::endl;
964 for (
size_t i = 0; i < numLocalRows; i++) {
966 for (
size_t ll = D0rowptr[i]; ll < D0rowptr[i+1]; ll++) {
971 for (
size_t jj = Prowptr[l]; jj < Prowptr[l+1]; jj++) {
974 for (
size_t k = 0; k < dim; k++) {
976 SC n = nullspace[k][i];
978 if (P11_status[jNew] == ST_INVALID || P11_status[jNew] < nnz_old) {
979 P11_status[jNew] = nnz;
980 P11colind[nnz] = jNew;
981 P11vals[nnz] = 0.5 * v * n;
986 P11vals[P11_status[jNew]] += 0.5 * v * n;
995 P11rowptr[numLocalRows] = nnz;
999 for (
size_t i = 0; i < numLocalRows; i++) {
1001 for (
size_t ll = D0rowptr[i]; ll < D0rowptr[i+1]; ll++) {
1002 LO l = D0colind[ll];
1003 for (
size_t jj = Prowptr[l]; jj < Prowptr[l+1]; jj++) {
1006 for (
size_t k = 0; k < dim; k++) {
1008 SC n = nullspace[k][i];
1010 if (P11_status[jNew] == ST_INVALID || P11_status[jNew] < nnz_old) {
1011 P11_status[jNew] = nnz;
1012 P11colind[nnz] = jNew;
1013 P11vals[nnz] = 0.5 * v * n;
1018 P11vals[P11_status[jNew]] += 0.5 * v * n;
1027 P11rowptr[numLocalRows] = nnz;
1035 P11colind_RCP.
resize(nnz);
1038 P11Crs->setAllValues(P11rowptr_RCP, P11colind_RCP, P11vals_RCP);
1039 P11Crs->expertStaticFillComplete(blockDomainMap, SM_Matrix_->getRangeMap());
1041 TEUCHOS_TEST_FOR_EXCEPTION(
false,std::invalid_argument,algo <<
" is not a valid option for \"refmaxwell: prolongator compute algorithm\"");
1042 #ifdef HAVE_MUELU_KOKKOS_REFACTOR 1047 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1053 if (parameterList_.get<
bool>(
"rap: triple product",
false) ==
false) {
1061 MultiplyRAP(*P11_,
true, *SM_Matrix_,
false, *P11_,
false, *Matrix1,
true,
true);
1063 if (parameterList_.get<
bool>(
"rap: fix zero diagonals",
true))
1066 if(disable_addon_==
true) {
1074 "MueLu::RefMaxwell::formCoarseMatrix(): Inverse of " 1075 "lumped mass matrix required for add-on (i.e. M0inv_Matrix is null)");
1088 if (M0inv_Matrix_->getGlobalMaxNumRowEntries()<=1) {
1091 #ifdef HAVE_MUELU_KOKKOS_REFACTOR 1094 ZT = Utilities_kokkos::Transpose(*Z);
1100 RCP<Vector> diag = VectorFactory::Build(M0inv_Matrix_->getRowMap());
1101 M0inv_Matrix_->getLocalDiagCopy(*diag);
1102 if (Z->getRowMap()->isSameAs(*(diag->getMap())))
1103 Z->leftScale(*diag);
1105 RCP<Import> importer = ImportFactory::Build(diag->getMap(),Z->getRowMap());
1106 RCP<Vector> diag2 = VectorFactory::Build(Z->getRowMap());
1108 Z->leftScale(*diag2);
1111 }
else if (parameterList_.get<
bool>(
"rap: triple product",
false) ==
false) {
1120 MultiplyRAP(*Z,
true, *M0inv_Matrix_,
false, *Z,
false, *Matrix2,
true,
true);
1123 Xpetra::MatrixMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::TwoMatrixAdd(*Matrix1,
false,(Scalar)1.0,*Matrix2,
false,(Scalar)1.0,AH_,GetOStream(
Runtime0));
1124 AH_->fillComplete();
1128 size_t dim = Nullspace_->getNumVectors();
1129 AH_->SetFixedBlockSize(dim);
1130 AH_->setObjectLabel(
"RefMaxwell (1,1)");
1135 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1137 SM_Matrix_ = SM_Matrix_new;
1141 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1147 residual_->update(one, RHS, negone);
1152 HierarchyH_->Iterate(*P11res_, *P11x_, 1,
true);
1153 Hierarchy22_->Iterate(*D0res_, *D0x_, 1,
true);
1158 X.
update(one, *residual_, one);
1163 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1169 residual_->update(one, RHS, negone);
1171 HierarchyH_->Iterate(*P11res_, *P11x_, 1,
true);
1173 X.
update(one, *residual_, one);
1177 residual_->update(one, RHS, negone);
1179 Hierarchy22_->Iterate(*D0res_, *D0x_, 1,
true);
1181 X.
update(one, *residual_, one);
1185 residual_->update(one, RHS, negone);
1187 HierarchyH_->Iterate(*P11res_, *P11x_, 1,
true);
1189 X.
update(one, *residual_, one);
1194 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1200 residual_->update(one, RHS, negone);
1202 Hierarchy22_->Iterate(*D0res_, *D0x_, 1,
true);
1204 X.
update(one, *residual_, one);
1208 residual_->update(one, RHS, negone);
1210 HierarchyH_->Iterate(*P11res_, *P11x_, 1,
true);
1212 X.
update(one, *residual_, one);
1216 residual_->update(one, RHS, negone);
1218 Hierarchy22_->Iterate(*D0res_, *D0x_, 1,
true);
1220 X.
update(one, *residual_, one);
1224 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1230 residual_->update(one, RHS, negone);
1234 HierarchyH_->Iterate(*P11res_, *P11x_, 1,
true);
1238 X.
update(one, *residual_, one);
1243 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1247 Scalar beta)
const {
1253 P11res_ = MultiVectorFactory::Build(P11_->getDomainMap(),X.
getNumVectors());
1254 P11x_ = MultiVectorFactory::Build(P11_->getDomainMap(),X.
getNumVectors());
1255 D0res_ = MultiVectorFactory::Build(D0_Matrix_->getDomainMap(),X.
getNumVectors());
1256 D0x_ = MultiVectorFactory::Build(D0_Matrix_->getDomainMap(),X.
getNumVectors());
1257 residual_ = MultiVectorFactory::Build(SM_Matrix_->getDomainMap(),1);
1261 #if defined(HAVE_MUELU_IFPACK2) && (!defined(HAVE_MUELU_EPETRA) || defined(HAVE_MUELU_INST_DOUBLE_INT_INT)) 1262 if (useHiptmairSmoothing_) {
1265 hiptmairPreSmoother_->apply(tRHS, tX);
1269 PreSmoother_->Apply(X, RHS, use_as_preconditioner_);
1272 if(mode_==
"additive")
1273 applyInverseAdditive(RHS,X);
1274 else if(mode_==
"121")
1275 applyInverse121(RHS,X);
1276 else if(mode_==
"212")
1277 applyInverse212(RHS,X);
1278 else if(mode_==
"11only")
1279 applyInverse11only(RHS,X);
1280 else if(mode_==
"none") {
1284 applyInverseAdditive(RHS,X);
1287 #if defined(HAVE_MUELU_IFPACK2) && (!defined(HAVE_MUELU_EPETRA) || defined(HAVE_MUELU_INST_DOUBLE_INT_INT)) 1288 if (useHiptmairSmoothing_)
1292 hiptmairPostSmoother_->apply(tRHS, tX);
1296 PostSmoother_->Apply(X, RHS,
false);
1301 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1307 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1320 HierarchyH_ = Teuchos::null;
1321 Hierarchy22_ = Teuchos::null;
1322 PreSmoother_ = Teuchos::null;
1323 PostSmoother_ = Teuchos::null;
1324 disable_addon_ =
false;
1328 setParameters(List);
1330 D0_Matrix_ = D0_Matrix;
1331 M0inv_Matrix_ = M0inv_Matrix;
1332 M1_Matrix_ = M1_Matrix;
1334 Nullspace_ = Nullspace;
1337 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1340 out <<
"\n--------------------------------------------------------------------------------\n" <<
1341 "--- RefMaxwell Summary ---\n" 1342 "--------------------------------------------------------------------------------" << std::endl;
1345 GlobalOrdinal numRows;
1348 numRows = SM_Matrix_->getGlobalNumRows();
1349 nnz = SM_Matrix_->getGlobalNumEntries();
1352 int rowspacer = 3;
while (tt != 0) { tt /= 10; rowspacer++; }
1354 int nnzspacer = 2;
while (tt != 0) { tt /= 10; nnzspacer++; }
1356 out <<
"block " << std::setw(rowspacer) <<
" rows " << std::setw(nnzspacer) <<
" nnz " << std::setw(9) <<
" nnz/row" << std::endl;
1357 out <<
"(1, 1)" << std::setw(rowspacer) << numRows << std::setw(nnzspacer) << nnz << std::setw(9) << as<double>(nnz) / numRows << std::endl;
1359 numRows = A22_->getGlobalNumRows();
1360 nnz = A22_->getGlobalNumEntries();
1362 out <<
"(2, 2)" << std::setw(rowspacer) << numRows << std::setw(nnzspacer) << nnz << std::setw(9) << as<double>(nnz) / numRows << std::endl;
1370 #define MUELU_REFMAXWELL_SHORT 1371 #endif //ifdef MUELU_REFMAXWELL_DEF_HPP Important warning messages (one line)
Generic Smoother Factory for generating the smoothers of the MG hierarchy.
This class specifies the default factory that should generate some data on a Level if the data does n...
static void CheckRepairMainDiagonal(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >> &Ac, bool const &repairZeroDiagonals, Teuchos::FancyOStream &fos)
Various adapters that will create a MueLu preconditioner that is an Xpetra::Matrix.
void parallel_for(const ExecPolicy &policy, const FunctorType &functor, const std::string &str="", typename Impl::enable_if< Kokkos::Impl::is_execution_policy< ExecPolicy >::value >::type *=0)
T & Get(const std::string &ename, const FactoryBase *factory=NoFactory::get())
Get data without decrementing associated storage counter (i.e., read-only access). Usage: Level->Get< RCP<Matrix> >("A", factory) if factory == NULL => use default factory.
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.
Factory for generating coarse level map. Used by TentativePFactory.
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 magnitudeType eps()
T & get(const std::string &name, T def_value)
Class that encapsulates external library smoothers.
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
static void MultiplyRAP(const Matrix &R, bool transposeR, const Matrix &A, bool transposeA, const Matrix &P, bool transposeP, Matrix &Ac, bool call_FillComplete_on_result=true, bool doOptimizeStorage=true, const std::string &label=std::string(), const RCP< ParameterList > ¶ms=null)
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)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
static void Write(const std::string &fileName, const Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > &M)
bool isSublist(const std::string &name) const
void ZeroDirichletCols(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Kokkos::View< const bool *, typename Node::device_type > &dirichletCols, Scalar replaceWith)
Teuchos::RCP< MueLu::Hierarchy< Scalar, LocalOrdinal, GlobalOrdinal, Node > > CreateXpetraPreconditioner(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > op, const Teuchos::ParameterList &inParamList, Teuchos::RCP< Xpetra::MultiVector< double, LocalOrdinal, GlobalOrdinal, Node > > coords=Teuchos::null, Teuchos::RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > nullspace=Teuchos::null)
Helper function to create a MueLu preconditioner that can be used by Xpetra.Given an Xpetra::Matrix...
void ZeroDirichletRows(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Kokkos::View< const bool *, typename Node::device_type > &dirichletRows, Scalar replaceWith)
One-liner description of what is happening.
Namespace for MueLu classes and methods.
void SetPreviousLevel(const RCP< Level > &previousLevel)
void SetFactoryManager(const RCP< const FactoryManagerBase > &factoryManager)
Set default factories (used internally by Hierarchy::SetLevel()).
void setParameters(Teuchos::ParameterList &list)
Set parameters.
static RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > RealValuedToScalarMultiVector(RCP< Xpetra::MultiVector< double, LocalOrdinal, GlobalOrdinal, Node > > X)
void applyInverseAdditive(const MultiVector &RHS, MultiVector &X) const
apply additive algorithm for 2x2 solve
void buildProlongator()
Setup the prolongator for the (1,1)-block.
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec)
Teuchos::ScalarTraits< Scalar >::magnitudeType magnitudeType
Factory for building tentative prolongator.
virtual void update(const Scalar &alpha, const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, const Scalar &beta)=0
void setlib(Xpetra::UnderlyingLib lib2)
static RCP< Time > getNewTimer(const std::string &name)
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
void resize(const size_type n, const T &val=T())
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::VERB_HIGH) const
void applyInverse121(const MultiVector &RHS, MultiVector &X) const
apply 1-2-1 algorithm for 2x2 solve
void Build(Level ¤tLevel) const
Build an object with this factory.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
virtual void SetParameterList(const ParameterList ¶mList)
Set parameters from a parameter list and return with default values.
Class that holds all level-specific information.
void compute()
Setup the preconditioner.
Teuchos::RCP< const Map > getRangeMap() const
Returns the Xpetra::Map object associated with the range of this operator.
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 > ¶ms=Teuchos::null)
virtual void setObjectLabel(const std::string &objectLabel)
AmalgamationFactory for subblocks of strided map based amalgamation data.
void initialize(const Teuchos::RCP< Matrix > &D0_Matrix, const Teuchos::RCP< Matrix > &M0inv_Matrix, const Teuchos::RCP< Matrix > &M1_Matrix, const Teuchos::RCP< MultiVector > &Nullspace, const Teuchos::RCP< RealValuedMultiVector > &Coords, Teuchos::ParameterList &List)
static void ZeroDirichletRows(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const std::vector< LocalOrdinal > &dirichletRows, Scalar replaceWith=Teuchos::ScalarTraits< Scalar >::zero())
void Set(const std::string &ename, const T &entry, const FactoryBase *factory=NoFactory::get())
static void TwoMatrixAdd(const Matrix &A, bool transposeA, SC alpha, Matrix &B, SC beta)
void SetParameter(const std::string &name, const ParameterEntry &entry)
Set a parameter directly as a ParameterEntry.
bool hasTransposeApply() const
Indicates whether this operator supports applying the adjoint operator.
void Build(Level ¤tLevel) const
Creates pre and post smoothers.
void applyInverse212(const MultiVector &RHS, MultiVector &X) const
apply 2-1-2 algorithm for 2x2 solve
Factory for creating a graph base on a given matrix.
static void Multiply(const Matrix &A, bool transposeA, const Matrix &B, bool transposeB, Matrix &C, bool call_FillComplete_on_result=true, bool doOptimizeStorage=true, const std::string &label=std::string(), const RCP< ParameterList > ¶ms=null)
void SetLevelID(int levelID)
Set level number.
void deep_copy(const View< DT, DP... > &dst, typename ViewTraits< DT, DP... >::const_value_type &value, typename std::enable_if< std::is_same< typename ViewTraits< DT, DP... >::specialize, void >::value >::type *=0)
static Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Read(const std::string &fileName, Xpetra::UnderlyingLib lib, const RCP< const Teuchos::Comm< int > > &comm, bool binary=false)
Class for transferring coordinates from a finer level to a coarser one.
Factory for building restriction operators.
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)
ParameterList & sublist(const std::string &name, bool mustAlreadyExist=false, const std::string &docString="")
virtual void SetFactory(const std::string &varName, const RCP< const FactoryBase > &factory)
Configuration.
void formCoarseMatrix()
Compute P11^{T}*A*P11 efficiently.
void applyInverse11only(const MultiVector &RHS, MultiVector &X) const
apply solve to 1-1 block only
void apply(const MultiVector &X, MultiVector &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=Teuchos::ScalarTraits< Scalar >::one(), Scalar beta=Teuchos::ScalarTraits< Scalar >::zero()) const
#define TEUCHOS_ASSERT(assertion_test)
static void ZeroDirichletCols(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Teuchos::ArrayRCP< const bool > &dirichletCols, Scalar replaceWith=Teuchos::ScalarTraits< Scalar >::zero())
#define TEUCHOS_ASSERT_EQUALITY(val1, val2)
Factory for building coarse matrices.
virtual size_t getNumVectors() const=0
static RCP< Tpetra::RowMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraRow(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
Factory for building uncoupled aggregates.
static Teuchos::RCP< Map< LocalOrdinal, GlobalOrdinal, Node > > Build(UnderlyingLib lib, global_size_t numGlobalElements, GlobalOrdinal indexBase, const Teuchos::RCP< const Teuchos::Comm< int > > &comm, LocalGlobal lg=Xpetra::GloballyDistributed, const Teuchos::RCP< Node > &=Teuchos::null)
Factory for building a thresholded operator.
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)
void resetMatrix(Teuchos::RCP< Matrix > SM_Matrix_new)
Reset system matrix.
Teuchos::RCP< const Map > getDomainMap() const
Returns the Xpetra::Map object associated with the domain of this operator.
void Request(const FactoryBase &factory)
Increment the storage counter for all the inputs of a factory.