46 #ifndef MUELU_PGPFACTORY_DEF_HPP 47 #define MUELU_PGPFACTORY_DEF_HPP 51 #include <Xpetra_Vector.hpp> 52 #include <Xpetra_MultiVectorFactory.hpp> 53 #include <Xpetra_VectorFactory.hpp> 54 #include <Xpetra_Import.hpp> 55 #include <Xpetra_ImportFactory.hpp> 56 #include <Xpetra_Export.hpp> 57 #include <Xpetra_ExportFactory.hpp> 58 #include <Xpetra_Matrix.hpp> 59 #include <Xpetra_MatrixMatrix.hpp> 63 #include "MueLu_PerfUtils.hpp" 66 #include "MueLu_SmootherFactory.hpp" 67 #include "MueLu_TentativePFactory.hpp" 68 #include "MueLu_Utilities.hpp" 72 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
74 RCP<ParameterList> validParamList = rcp(
new ParameterList());
76 validParamList->set< RCP<const FactoryBase> >(
"A", Teuchos::null,
"Generating factory of the matrix A used during the prolongator smoothing process");
77 validParamList->set< RCP<const FactoryBase> >(
"P", Teuchos::null,
"Tentative prolongator factory");
79 validParamList->set<
bool > (
"ReUseRowBasedOmegas",
false,
"Reuse omegas for prolongator for restrictor");
81 return validParamList;
84 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
86 SetParameter(
"Minimization norm", ParameterEntry(minnorm));
89 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
91 const ParameterList& pL = GetParameterList();
95 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
97 Input(fineLevel,
"A");
102 RCP<const FactoryBase> initialPFact = GetFactory(
"P");
103 if (initialPFact == Teuchos::null) { initialPFact = coarseLevel.
GetFactoryManager()->GetFactory(
"Ptent"); }
104 coarseLevel.
DeclareInput(
"P", initialPFact.get(),
this);
120 const ParameterList & pL = GetParameterList();
121 bool bReUseRowBasedOmegas = pL.get<
bool>(
"ReUseRowBasedOmegas");
122 if( bReUseRowBasedOmegas ==
true && restrictionMode_ ==
true ) {
127 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
129 FactoryMonitor m(*
this,
"Prolongator smoothing (PG-AMG)", coarseLevel);
132 RCP<Matrix> A = Get< RCP<Matrix> >(fineLevel,
"A");
137 RCP<const FactoryBase> initialPFact = GetFactory(
"P");
138 if (initialPFact == Teuchos::null) { initialPFact = coarseLevel.
GetFactoryManager()->GetFactory(
"Ptent"); }
139 RCP<Matrix> Ptent = coarseLevel.
Get< RCP<Matrix> >(
"P", initialPFact.get());
142 if(restrictionMode_) {
148 bool doFillComplete=
true;
149 bool optimizeStorage=
true;
150 RCP<Matrix> DinvAP0 = Xpetra::MatrixMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Multiply(*A,
false, *Ptent,
false, GetOStream(
Statistics2), doFillComplete, optimizeStorage);
153 optimizeStorage=
false;
159 Teuchos::RCP<Vector > RowBasedOmega = Teuchos::null;
161 const ParameterList & pL = GetParameterList();
162 bool bReUseRowBasedOmegas = pL.get<
bool>(
"ReUseRowBasedOmegas");
163 if(restrictionMode_ ==
false || bReUseRowBasedOmegas ==
false) {
166 ComputeRowBasedOmega(fineLevel, coarseLevel, A, Ptent, DinvAP0, RowBasedOmega);
170 RowBasedOmega = coarseLevel.
Get<Teuchos::RCP<Vector > >(
"RowBasedOmega",
this);
178 Teuchos::RCP<const Export> exporter =
179 ExportFactory::Build(RowBasedOmega->getMap(), A->getRangeMap());
181 Teuchos::RCP<Vector > noRowBasedOmega =
182 VectorFactory::Build(A->getRangeMap());
184 noRowBasedOmega->doExport(*RowBasedOmega, *exporter, Xpetra::INSERT);
189 Teuchos::RCP<const Import > importer =
190 ImportFactory::Build(A->getRangeMap(), A->getRowMap());
193 RowBasedOmega->doImport(*noRowBasedOmega, *importer, Xpetra::INSERT);
196 Teuchos::ArrayRCP< Scalar > RowBasedOmega_local = RowBasedOmega->getDataNonConst(0);
199 RCP<Matrix> P_smoothed = Teuchos::null;
202 Xpetra::MatrixMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::TwoMatrixAdd(*Ptent,
false, Teuchos::ScalarTraits<Scalar>::one(),
203 *DinvAP0,
false, -Teuchos::ScalarTraits<Scalar>::one(),
205 P_smoothed->fillComplete(Ptent->getDomainMap(), Ptent->getRangeMap());
209 RCP<ParameterList> params = rcp(
new ParameterList());
210 params->set(
"printLoadBalancingInfo",
true);
213 if (!restrictionMode_) {
215 Set(coarseLevel,
"P", P_smoothed);
221 if (Ptent->IsView(
"stridedMaps"))
222 P_smoothed->CreateView(
"stridedMaps", Ptent);
227 Set(coarseLevel,
"R", R);
233 if (Ptent->IsView(
"stridedMaps"))
234 R->CreateView(
"stridedMaps", Ptent,
true);
239 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
241 FactoryMonitor m(*
this,
"PgPFactory::ComputeRowBasedOmega", coarseLevel);
243 typedef typename Teuchos::ScalarTraits<Scalar>::magnitudeType Magnitude;
244 Scalar sZero = Teuchos::ScalarTraits<Scalar>::zero();
245 Magnitude mZero = Teuchos::ScalarTraits<Scalar>::magnitude(sZero);
247 Teuchos::RCP<Vector > Numerator = Teuchos::null;
248 Teuchos::RCP<Vector > Denominator = Teuchos::null;
250 const ParameterList & pL = GetParameterList();
267 bool doFillComplete=
true;
268 bool optimizeStorage=
false;
269 RCP<Matrix> AP0 = Xpetra::MatrixMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Multiply(*A,
false, *P0,
false, GetOStream(
Statistics2), doFillComplete, optimizeStorage);
272 RCP<Matrix> ADinvAP0 = Xpetra::MatrixMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Multiply(*A,
false, *DinvAP0,
false, GetOStream(
Statistics2), doFillComplete, optimizeStorage);
274 Numerator = VectorFactory::Build(ADinvAP0->getColMap(),
true);
275 Denominator = VectorFactory::Build(ADinvAP0->getColMap(),
true);
276 MultiplyAll(AP0, ADinvAP0, Numerator);
277 MultiplySelfAll(ADinvAP0, Denominator);
288 Numerator = VectorFactory::Build(DinvAP0->getColMap(),
true);
289 Denominator = VectorFactory::Build(DinvAP0->getColMap(),
true);
290 MultiplyAll(P0, DinvAP0, Numerator);
291 MultiplySelfAll(DinvAP0, Denominator);
305 bool doFillComplete=
true;
306 bool optimizeStorage=
true;
308 RCP<Matrix> DinvADinvAP0 = Xpetra::MatrixMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Multiply(*A,
false, *DinvAP0,
false, GetOStream(
Statistics2), doFillComplete, optimizeStorage);
310 diagA = Teuchos::ArrayRCP<Scalar>();
312 Numerator = VectorFactory::Build(DinvADinvAP0->getColMap(),
true);
313 Denominator = VectorFactory::Build(DinvADinvAP0->getColMap(),
true);
314 MultiplyAll(DinvAP0, DinvADinvAP0, Numerator);
315 MultiplySelfAll(DinvADinvAP0, Denominator);
319 TEUCHOS_TEST_FOR_EXCEPTION(
true,
Exceptions::RuntimeError,
"MueLu::PgPFactory::Build: minimization mode not supported. error");
324 Teuchos::RCP<Vector > ColBasedOmega =
325 VectorFactory::Build(Numerator->getMap(),
true);
327 ColBasedOmega->putScalar(-666);
329 Teuchos::ArrayRCP< const Scalar > Numerator_local = Numerator->getData(0);
330 Teuchos::ArrayRCP< const Scalar > Denominator_local = Denominator->getData(0);
331 Teuchos::ArrayRCP< Scalar > ColBasedOmega_local = ColBasedOmega->getDataNonConst(0);
334 Magnitude min_local = 1000000.0;
335 Magnitude max_local = 0.0;
336 for(
LocalOrdinal i = 0; i < Teuchos::as<LocalOrdinal>(Numerator->getLocalLength()); i++) {
337 if(Teuchos::ScalarTraits<Scalar>::magnitude(Denominator_local[i]) == mZero)
339 ColBasedOmega_local[i] = 0.0;
344 ColBasedOmega_local[i] = Numerator_local[i] / Denominator_local[i];
347 if(Teuchos::ScalarTraits<Scalar>::magnitude(ColBasedOmega_local[i]) < mZero) {
348 ColBasedOmega_local[i] = Teuchos::ScalarTraits<Scalar>::zero();
356 if(Teuchos::ScalarTraits<Scalar>::magnitude(ColBasedOmega_local[i]) >= 0.8) {
357 ColBasedOmega_local[i] = 0.0;
360 if(Teuchos::ScalarTraits<Scalar>::magnitude(ColBasedOmega_local[i]) < min_local)
361 { min_local = Teuchos::ScalarTraits<Scalar>::magnitude(ColBasedOmega_local[i]); }
362 if(Teuchos::ScalarTraits<Scalar>::magnitude(ColBasedOmega_local[i]) > max_local)
363 { max_local = Teuchos::ScalarTraits<Scalar>::magnitude(ColBasedOmega_local[i]); }
371 MueLu_sumAll(A->getRowMap()->getComm(), zero_local, zero_all);
372 MueLu_sumAll(A->getRowMap()->getComm(), nan_local, nan_all);
373 MueLu_minAll(A->getRowMap()->getComm(), min_local, min_all);
374 MueLu_maxAll(A->getRowMap()->getComm(), max_local, max_all);
381 default: GetOStream(
Statistics1) <<
"unknown)" << std::endl;
break;
383 GetOStream(
Statistics1) <<
"Damping parameter: min = " << min_all <<
", max = " << max_all << std::endl;
384 GetOStream(
Statistics) <<
"# negative omegas: " << zero_all <<
" out of " << ColBasedOmega->getGlobalLength() <<
" column-based omegas" << std::endl;
385 GetOStream(
Statistics) <<
"# NaNs: " << nan_all <<
" out of " << ColBasedOmega->getGlobalLength() <<
" column-based omegas" << std::endl;
388 if(coarseLevel.
IsRequested(
"ColBasedOmega",
this)) {
389 coarseLevel.
Set(
"ColBasedOmega", ColBasedOmega,
this);
395 VectorFactory::Build(DinvAP0->getRowMap(),
true);
397 RowBasedOmega->putScalar(-666);
399 bool bAtLeastOneDefined =
false;
400 Teuchos::ArrayRCP< Scalar > RowBasedOmega_local = RowBasedOmega->getDataNonConst(0);
401 for(
LocalOrdinal row = 0; row<Teuchos::as<LocalOrdinal>(A->getNodeNumRows()); row++) {
402 Teuchos::ArrayView<const LocalOrdinal> lindices;
403 Teuchos::ArrayView<const Scalar> lvals;
404 DinvAP0->getLocalRowView(row, lindices, lvals);
405 bAtLeastOneDefined =
false;
406 for(
size_t j=0; j<Teuchos::as<size_t>(lindices.size()); j++) {
407 Scalar omega = ColBasedOmega_local[lindices[j]];
408 if (Teuchos::ScalarTraits<Scalar>::magnitude(omega) != -666) {
409 bAtLeastOneDefined =
true;
410 if(Teuchos::ScalarTraits<Scalar>::magnitude(RowBasedOmega_local[row]) == -666) RowBasedOmega_local[row] = omega;
411 else if(Teuchos::ScalarTraits<Scalar>::magnitude(omega) < Teuchos::ScalarTraits<Scalar>::magnitude(RowBasedOmega_local[row])) RowBasedOmega_local[row] = omega;
414 if(bAtLeastOneDefined ==
true) {
415 if(Teuchos::ScalarTraits<Scalar>::magnitude(RowBasedOmega_local[row]) < mZero)
416 RowBasedOmega_local[row] = sZero;
420 if(coarseLevel.
IsRequested(
"RowBasedOmega",
this)) {
421 Set(coarseLevel,
"RowBasedOmega", RowBasedOmega);
425 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
429 TEUCHOS_TEST_FOR_EXCEPTION(!InnerProdVec->getMap()->isSameAs(*Op->getColMap()),
Exceptions::RuntimeError,
"MueLu::PgPFactory::MultiplySelfAll: map of InnerProdVec must be same as column map of operator. error");
431 Teuchos::ArrayRCP< Scalar > InnerProd_local = InnerProdVec->getDataNonConst(0);
433 Teuchos::ArrayView<const LocalOrdinal> lindices;
434 Teuchos::ArrayView<const Scalar> lvals;
436 for(
size_t n=0; n<Op->getNodeNumRows(); n++) {
437 Op->getLocalRowView(n, lindices, lvals);
438 for(
size_t i=0; i<Teuchos::as<size_t>(lindices.size()); i++) {
439 InnerProd_local[lindices[i]] += lvals[i]*lvals[i];
442 InnerProd_local = Teuchos::ArrayRCP< Scalar >();
445 Teuchos::RCP<const Export> exporter =
446 ExportFactory::Build(Op->getColMap(), Op->getDomainMap());
448 Teuchos::RCP<Vector > nonoverlap =
449 VectorFactory::Build(Op->getDomainMap());
451 nonoverlap->doExport(*InnerProdVec, *exporter, Xpetra::ADD);
456 Teuchos::RCP<const Import > importer =
457 ImportFactory::Build(Op->getDomainMap(), Op->getColMap());
460 InnerProdVec->doImport(*nonoverlap, *importer, Xpetra::INSERT);
464 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
467 TEUCHOS_TEST_FOR_EXCEPTION(!left->getDomainMap()->isSameAs(*right->getDomainMap()),
Exceptions::RuntimeError,
"MueLu::PgPFactory::MultiplyAll: domain maps of left and right do not match. Error.");
468 TEUCHOS_TEST_FOR_EXCEPTION(!left->getRowMap()->isSameAs(*right->getRowMap()),
Exceptions::RuntimeError,
"MueLu::PgPFactory::MultiplyAll: row maps of left and right do not match. Error.");
469 #if 1 // 1=new "fast code, 0=old "slow", but safe code 470 #if 0 // not necessary - remove me 471 if(InnerProdVec->getMap()->isSameAs(*left->getColMap())) {
474 std::vector<LocalOrdinal> NewRightLocal(right->getColMap()->getNodeNumElements(), Teuchos::as<LocalOrdinal>(left->getColMap()->getNodeNumElements()+1));
477 for (
size_t j=0; j < right->getColMap()->getNodeNumElements(); j++) {
478 while ( (i < Teuchos::as<LocalOrdinal>(left->getColMap()->getNodeNumElements())) &&
479 (left->getColMap()->getGlobalElement(i) < right->getColMap()->getGlobalElement(j)) ) i++;
480 if (left->getColMap()->getGlobalElement(i) == right->getColMap()->getGlobalElement(j)) {
481 NewRightLocal[j] = i;
485 Teuchos::ArrayRCP< Scalar > InnerProd_local = InnerProdVec->getDataNonConst(0);
486 std::vector<Scalar> temp_array(left->getColMap()->getNodeNumElements()+1, 0.0);
488 for(
size_t n=0; n<right->getNodeNumRows(); n++) {
489 Teuchos::ArrayView<const LocalOrdinal> lindices_left;
490 Teuchos::ArrayView<const Scalar> lvals_left;
491 Teuchos::ArrayView<const LocalOrdinal> lindices_right;
492 Teuchos::ArrayView<const Scalar> lvals_right;
494 left->getLocalRowView (n, lindices_left, lvals_left);
495 right->getLocalRowView(n, lindices_right, lvals_right);
497 for(
size_t j=0; j<Teuchos::as<size_t>(lindices_right.size()); j++) {
498 temp_array[NewRightLocal[lindices_right[j] ] ] = lvals_right[j];
500 for (
size_t j=0; j < Teuchos::as<size_t>(lindices_left.size()); j++) {
501 InnerProd_local[lindices_left[j]] += temp_array[lindices_left[j] ]*lvals_left[j];
503 for (
size_t j=0; j < Teuchos::as<size_t>(lindices_right.size()); j++) {
504 temp_array[NewRightLocal[lindices_right[j] ] ] = 0.0;
508 Teuchos::RCP<const Export> exporter =
509 ExportFactory::Build(left->getColMap(), left->getDomainMap());
511 Teuchos::RCP<Vector > nonoverlap =
512 VectorFactory::Build(left->getDomainMap());
514 nonoverlap->doExport(*InnerProdVec, *exporter, Xpetra::ADD);
519 Teuchos::RCP<const Import > importer =
520 ImportFactory::Build(left->getDomainMap(), left->getColMap());
523 InnerProdVec->doImport(*nonoverlap, *importer, Xpetra::INSERT);
527 #endif // end remove me 528 if(InnerProdVec->getMap()->isSameAs(*right->getColMap())) {
529 size_t szNewLeftLocal = TEUCHOS_MAX(left->getColMap()->getNodeNumElements(), right->getColMap()->getNodeNumElements());
530 Teuchos::RCP<std::vector<LocalOrdinal> > NewLeftLocal = Teuchos::rcp(
new std::vector<LocalOrdinal>(szNewLeftLocal, Teuchos::as<LocalOrdinal>(right->getColMap()->getMaxLocalIndex()+1)));
533 for (
size_t i=0; i < left->getColMap()->getNodeNumElements(); i++) {
534 while ( (j < Teuchos::as<LocalOrdinal>(right->getColMap()->getNodeNumElements())) &&
535 (right->getColMap()->getGlobalElement(j) < left->getColMap()->getGlobalElement(i)) ) j++;
536 if (right->getColMap()->getGlobalElement(j) == left->getColMap()->getGlobalElement(i)) {
537 (*NewLeftLocal)[i] = j;
545 Teuchos::ArrayRCP< Scalar > InnerProd_local = InnerProdVec->getDataNonConst(0);
546 Teuchos::RCP<std::vector<Scalar> > temp_array = Teuchos::rcp(
new std::vector<Scalar>(right->getColMap()->getMaxLocalIndex()+2, 0.0));
548 for(
size_t n=0; n<left->getNodeNumRows(); n++) {
549 Teuchos::ArrayView<const LocalOrdinal> lindices_left;
550 Teuchos::ArrayView<const Scalar> lvals_left;
551 Teuchos::ArrayView<const LocalOrdinal> lindices_right;
552 Teuchos::ArrayView<const Scalar> lvals_right;
554 left->getLocalRowView (n, lindices_left, lvals_left);
555 right->getLocalRowView(n, lindices_right, lvals_right);
557 for(
size_t i=0; i<Teuchos::as<size_t>(lindices_left.size()); i++) {
558 (*temp_array)[(*NewLeftLocal)[lindices_left[i] ] ] = lvals_left[i];
560 for (
size_t i=0; i < Teuchos::as<size_t>(lindices_right.size()); i++) {
561 InnerProd_local[lindices_right[i]] += (*temp_array)[lindices_right[i] ] * lvals_right[i];
563 for (
size_t i=0; i < Teuchos::as<size_t>(lindices_left.size()); i++) {
564 (*temp_array)[(*NewLeftLocal)[lindices_left[i] ] ] = 0.0;
567 InnerProd_local = Teuchos::ArrayRCP< Scalar >();
569 Teuchos::RCP<const Export> exporter =
570 ExportFactory::Build(right->getColMap(), right->getDomainMap());
572 Teuchos::RCP<Vector> nonoverlap =
573 VectorFactory::Build(right->getDomainMap());
575 nonoverlap->doExport(*InnerProdVec, *exporter, Xpetra::ADD);
580 Teuchos::RCP<const Import > importer =
581 ImportFactory::Build(right->getDomainMap(), right->getColMap());
583 InnerProdVec->doImport(*nonoverlap, *importer, Xpetra::INSERT);
585 TEUCHOS_TEST_FOR_EXCEPTION(
true,
Exceptions::RuntimeError,
"MueLu::PgPFactory::MultiplyAll: map of InnerProdVec must be same as column map of left operator? Error.");
588 #else // old "safe" code 589 if(InnerProdVec->getMap()->isSameAs(*left->getColMap())) {
591 Teuchos::ArrayRCP< Scalar > InnerProd_local = InnerProdVec->getDataNonConst(0);
594 Teuchos::ArrayView<const LocalOrdinal> lindices_left;
595 Teuchos::ArrayView<const Scalar> lvals_left;
596 Teuchos::ArrayView<const LocalOrdinal> lindices_right;
597 Teuchos::ArrayView<const Scalar> lvals_right;
599 for(
size_t n=0; n<left->getNodeNumRows(); n++)
602 left->getLocalRowView (n, lindices_left, lvals_left);
603 right->getLocalRowView(n, lindices_right, lvals_right);
605 for(
size_t i=0; i<Teuchos::as<size_t>(lindices_left.size()); i++)
607 GlobalOrdinal left_gid = left->getColMap()->getGlobalElement(lindices_left[i]);
608 for(
size_t j=0; j<Teuchos::as<size_t>(lindices_right.size()); j++)
610 GlobalOrdinal right_gid= right->getColMap()->getGlobalElement(lindices_right[j]);
611 if(left_gid == right_gid)
613 InnerProd_local[lindices_left[i]] += lvals_left[i]*lvals_right[j];
621 Teuchos::RCP<const Export> exporter =
622 ExportFactory::Build(left->getColMap(), left->getDomainMap());
624 Teuchos::RCP<Vector > nonoverlap =
625 VectorFactory::Build(left->getDomainMap());
627 nonoverlap->doExport(*InnerProdVec, *exporter, Xpetra::ADD);
632 Teuchos::RCP<const Import > importer =
633 ImportFactory::Build(left->getDomainMap(), left->getColMap());
636 InnerProdVec->doImport(*nonoverlap, *importer, Xpetra::INSERT);
638 else if(InnerProdVec->getMap()->isSameAs(*right->getColMap())) {
639 Teuchos::ArrayRCP< Scalar > InnerProd_local = InnerProdVec->getDataNonConst(0);
641 Teuchos::ArrayView<const LocalOrdinal> lindices_left;
642 Teuchos::ArrayView<const Scalar> lvals_left;
643 Teuchos::ArrayView<const LocalOrdinal> lindices_right;
644 Teuchos::ArrayView<const Scalar> lvals_right;
646 for(
size_t n=0; n<left->getNodeNumRows(); n++)
648 left->getLocalRowView(n, lindices_left, lvals_left);
649 right->getLocalRowView(n, lindices_right, lvals_right);
651 for(
size_t i=0; i<Teuchos::as<size_t>(lindices_left.size()); i++)
653 GlobalOrdinal left_gid = left->getColMap()->getGlobalElement(lindices_left[i]);
654 for(
size_t j=0; j<Teuchos::as<size_t>(lindices_right.size()); j++)
656 GlobalOrdinal right_gid= right->getColMap()->getGlobalElement(lindices_right[j]);
657 if(left_gid == right_gid)
659 InnerProd_local[lindices_right[j]] += lvals_left[i]*lvals_right[j];
667 Teuchos::RCP<const Export> exporter =
668 ExportFactory::Build(right->getColMap(), right->getDomainMap());
670 Teuchos::RCP<Vector> nonoverlap =
671 VectorFactory::Build(right->getDomainMap());
673 nonoverlap->doExport(*InnerProdVec, *exporter, Xpetra::ADD);
678 Teuchos::RCP<const Import > importer =
679 ImportFactory::Build(right->getDomainMap(), right->getColMap());
682 InnerProdVec->doImport(*nonoverlap, *importer, Xpetra::INSERT);
685 TEUCHOS_TEST_FOR_EXCEPTION(
true,
Exceptions::RuntimeError,
"MueLu::PgPFactory::MultiplyAll: map of InnerProdVec must be same as column map of left or right operator? Error.");
690 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
692 std::cout <<
"TODO: remove me" << std::endl;
695 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
697 SetParameter(
"ReUseRowBasedOmegas", ParameterEntry(bReuse));
static Teuchos::ArrayRCP< Scalar > GetMatrixDiagonal(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.
#define MueLu_sumAll(rcpComm, in, out)
void ComputeRowBasedOmega(Level &fineLevel, Level &coarseLevel, const RCP< Matrix > &A, const RCP< Matrix > &P0, const RCP< Matrix > &DinvAP0, RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &RowBasedOmega) const
MueLu::DefaultLocalOrdinal LocalOrdinal
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 void MyOldScaleMatrix(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Teuchos::ArrayRCP< const Scalar > &scalingVector, bool doInverse=true, bool doFillComplete=true, bool doOptimizeStorage=true)
#define MueLu_maxAll(rcpComm, in, out)
Timer to be used in factories. Similar to Monitor but with additional timers.
Namespace for MueLu classes and methods.
#define MueLu_minAll(rcpComm, in, out)
Print even more statistics.
bool IsRequested(const std::string &ename, const FactoryBase *factory=NoFactory::get()) const
Test whether a need has been requested. Note: this tells nothing about whether the need's value exist...
void SetMinimizationMode(MinimizationNorm minnorm)
Set minimization mode (L2NORM for cheapest, ANORM more expensive, DINVANORM = default) ...
void Build(Level &fineLevel, Level &coarseLevel) const
Build method.
MueLu::DefaultScalar Scalar
MueLu::DefaultGlobalOrdinal GlobalOrdinal
Class that holds all level-specific information.
Timer to be used in factories. Similar to SubMonitor but adds a timer level by level.
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)
void MultiplySelfAll(const RCP< Matrix > &Op, Teuchos::RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &InnerProdVec) const
static std::string PrintMatrixInfo(const Matrix &A, const std::string &msgTag, RCP< const Teuchos::ParameterList > params=Teuchos::null)
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
MinimizationNorm GetMinimizationMode()
return minimization mode
void Set(const std::string &ename, const T &entry, const FactoryBase *factory=NoFactory::get())
void ReUseDampingParameters(bool bReuse)
void BuildP(Level &fineLevel, Level &coarseLevel) const
Abstract Build method.
const RCP< const FactoryManagerBase > GetFactoryManager()
returns the current factory manager
Exception throws to report errors in the internal logical of the program.
void DeclareInput(const std::string &ename, const FactoryBase *factory, const FactoryBase *requestedBy=NoFactory::get())
Callback from FactoryBase::CallDeclareInput() and FactoryBase::DeclareInput()
void MultiplyAll(const RCP< Matrix > &left, const RCP< Matrix > &right, Teuchos::RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &InnerProdVec) const