45 #ifndef MUELU_PARAMETERLISTINTERPRETER_DEF_HPP 46 #define MUELU_PARAMETERLISTINTERPRETER_DEF_HPP 48 #include <Teuchos_XMLParameterListHelpers.hpp> 50 #include <Xpetra_Matrix.hpp> 51 #include <Xpetra_MatrixUtils.hpp> 59 #include "MueLu_Hierarchy.hpp" 60 #include "MueLu_FactoryManager.hpp" 62 #include "MueLu_AggregationExportFactory.hpp" 63 #include "MueLu_AggregateQualityEstimateFactory.hpp" 64 #include "MueLu_BrickAggregationFactory.hpp" 65 #include "MueLu_CoalesceDropFactory.hpp" 66 #include "MueLu_CoarseMapFactory.hpp" 67 #include "MueLu_ConstraintFactory.hpp" 68 #include "MueLu_CoordinatesTransferFactory.hpp" 69 #include "MueLu_CoupledAggregationFactory.hpp" 70 #include "MueLu_DirectSolver.hpp" 71 #include "MueLu_EminPFactory.hpp" 73 #include "MueLu_FacadeClassFactory.hpp" 74 #include "MueLu_FactoryFactory.hpp" 75 #include "MueLu_FilteredAFactory.hpp" 76 #include "MueLu_GenericRFactory.hpp" 77 #include "MueLu_LineDetectionFactory.hpp" 79 #include "MueLu_NotayAggregationFactory.hpp" 80 #include "MueLu_NullspaceFactory.hpp" 81 #include "MueLu_PatternFactory.hpp" 82 #include "MueLu_PgPFactory.hpp" 83 #include "MueLu_RAPFactory.hpp" 84 #include "MueLu_RAPShiftFactory.hpp" 85 #include "MueLu_RebalanceAcFactory.hpp" 86 #include "MueLu_RebalanceTransferFactory.hpp" 87 #include "MueLu_RepartitionFactory.hpp" 88 #include "MueLu_SaPFactory.hpp" 89 #include "MueLu_ScaledNullspaceFactory.hpp" 90 #include "MueLu_SemiCoarsenPFactory.hpp" 91 #include "MueLu_SmootherFactory.hpp" 92 #include "MueLu_SmooVecCoalesceDropFactory.hpp" 93 #include "MueLu_TentativePFactory.hpp" 94 #include "MueLu_TogglePFactory.hpp" 95 #include "MueLu_ToggleCoordinatesTransferFactory.hpp" 96 #include "MueLu_TransPFactory.hpp" 97 #include "MueLu_UncoupledAggregationFactory.hpp" 98 #include "MueLu_HybridAggregationFactory.hpp" 99 #include "MueLu_ZoltanInterface.hpp" 100 #include "MueLu_Zoltan2Interface.hpp" 101 #include "MueLu_NodePartitionInterface.hpp" 103 #ifdef HAVE_MUELU_KOKKOS_REFACTOR 104 #include "MueLu_CoalesceDropFactory_kokkos.hpp" 105 #include "MueLu_CoarseMapFactory_kokkos.hpp" 106 #include "MueLu_CoordinatesTransferFactory_kokkos.hpp" 107 #include "MueLu_NullspaceFactory_kokkos.hpp" 108 #include "MueLu_SaPFactory_kokkos.hpp" 109 #include "MueLu_TentativePFactory_kokkos.hpp" 110 #include "MueLu_UncoupledAggregationFactory_kokkos.hpp" 113 #ifdef HAVE_MUELU_MATLAB 114 #include "../matlab/src/MueLu_MatlabSmoother_decl.hpp" 115 #include "../matlab/src/MueLu_MatlabSmoother_def.hpp" 116 #include "../matlab/src/MueLu_TwoLevelMatlabFactory_decl.hpp" 117 #include "../matlab/src/MueLu_TwoLevelMatlabFactory_def.hpp" 118 #include "../matlab/src/MueLu_SingleLevelMatlabFactory_decl.hpp" 119 #include "../matlab/src/MueLu_SingleLevelMatlabFactory_def.hpp" 122 #ifdef HAVE_MUELU_INTREPID2 123 #include "MueLu_IntrepidPCoarsenFactory.hpp" 126 #include <unordered_set> 130 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
132 RCP<Teuchos::TimeMonitor> tM = rcp(
new Teuchos::TimeMonitor(*Teuchos::TimeMonitor::getNewTimer(std::string(
"MueLu: ParameterListInterpreter (ParameterList)"))));
133 if(facadeFact == Teuchos::null)
134 facadeFact_ = Teuchos::rcp(
new FacadeClassFactory());
136 facadeFact_ = facadeFact;
138 if (paramList.isParameter(
"xml parameter file")) {
139 std::string filename = paramList.get(
"xml parameter file",
"");
140 if (filename.length() != 0) {
141 TEUCHOS_TEST_FOR_EXCEPTION(comm.is_null(), Exceptions::RuntimeError,
"xml parameter file requires a valid comm");
143 ParameterList paramList2 = paramList;
144 Teuchos::updateParametersFromXmlFileAndBroadcast(filename, Teuchos::Ptr<Teuchos::ParameterList>(¶mList2), *comm);
145 SetParameterList(paramList2);
148 SetParameterList(paramList);
152 SetParameterList(paramList);
156 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
158 RCP<Teuchos::TimeMonitor> tM = rcp(
new Teuchos::TimeMonitor(*Teuchos::TimeMonitor::getNewTimer(std::string(
"MueLu: ParameterListInterpreter (XML)"))));
159 if(facadeFact == Teuchos::null)
164 ParameterList paramList;
165 Teuchos::updateParametersFromXmlFileAndBroadcast(xmlFileName, Teuchos::Ptr<ParameterList>(¶mList), comm);
169 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
173 scalingFactor_= Teuchos::ScalarTraits<double>::one();
177 if (paramList.isSublist(
"Hierarchy")) {
178 SetFactoryParameterList(paramList);
180 }
else if (paramList.isParameter(
"MueLu preconditioner") ==
true) {
181 this->GetOStream(
Runtime0) <<
"Use facade class: " << paramList.get<std::string>(
"MueLu preconditioner") << std::endl;
182 Teuchos::RCP<ParameterList> pp = facadeFact_->SetParameterList(paramList);
183 SetFactoryParameterList(*pp);
187 ParameterList serialList, nonSerialList;
190 Validate(serialList);
191 SetEasyParameterList(paramList);
199 static inline bool areSame(
const ParameterList& list1,
const ParameterList& list2);
204 #define MUELU_SET_VAR_2LIST(paramList, defaultList, paramName, paramType, varName) \ 206 if (paramList.isParameter(paramName)) varName = paramList.get<paramType>(paramName); \ 207 else if (defaultList.isParameter(paramName)) varName = defaultList.get<paramType>(paramName); \ 208 else varName = MasterList::getDefault<paramType>(paramName); 210 #define MUELU_TEST_AND_SET_VAR(paramList, paramName, paramType, varName) \ 211 (paramList.isParameter(paramName) ? varName = paramList.get<paramType>(paramName), true : false) 215 #define MUELU_TEST_AND_SET_PARAM_2LIST(paramList, defaultList, paramName, paramType, listWrite) \ 217 if (paramList .isParameter(paramName)) listWrite.set(paramName, paramList .get<paramType>(paramName)); \ 218 else if (defaultList.isParameter(paramName)) listWrite.set(paramName, defaultList.get<paramType>(paramName)); \ 220 catch(Teuchos::Exceptions::InvalidParameterType&) { \ 221 TEUCHOS_TEST_FOR_EXCEPTION_PURE_MSG(true, Teuchos::Exceptions::InvalidParameterType, \ 222 "Error: parameter \"" << paramName << "\" must be of type " << Teuchos::TypeNameTraits<paramType>::name()); \ 225 #define MUELU_TEST_PARAM_2LIST(paramList, defaultList, paramName, paramType, cmpValue) \ 227 paramList.isParameter(paramName) ? paramList .get<paramType>(paramName) : ( \ 228 defaultList.isParameter(paramName) ? defaultList.get<paramType>(paramName) : \ 229 MasterList::getDefault<paramType>(paramName) ) ) ) 231 #ifndef HAVE_MUELU_KOKKOS_REFACTOR 232 #define MUELU_KOKKOS_FACTORY(varName, oldFactory, newFactory) \ 233 RCP<Factory> varName = rcp(new oldFactory()); 234 #define MUELU_KOKKOS_FACTORY_NO_DECL(varName, oldFactory, newFactory) \ 235 varName = rcp(new oldFactory()); 237 #define MUELU_KOKKOS_FACTORY(varName, oldFactory, newFactory) \ 238 RCP<Factory> varName; \ 239 if (!useKokkos_) varName = rcp(new oldFactory()); \ 240 else varName = rcp(new newFactory()); 241 #define MUELU_KOKKOS_FACTORY_NO_DECL(varName, oldFactory, newFactory) \ 242 if (!useKokkos_) varName = rcp(new oldFactory()); \ 243 else varName = rcp(new newFactory()); 246 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
249 ParameterList paramList;
251 MUELU_SET_VAR_2LIST(constParamList, constParamList,
"problem: type", std::string, problemType);
252 if (problemType !=
"unknown") {
254 paramList.setParameters(constParamList);
258 paramList = constParamList;
262 #if !defined(HAVE_MUELU_KOKKOS_REFACTOR) 265 # ifdef HAVE_MUELU_SERIAL 266 if (
typeid(
Node).name() ==
typeid(Kokkos::Compat::KokkosSerialWrapperNode).name())
269 # ifdef HAVE_MUELU_OPENMP 270 if (
typeid(
Node).name() ==
typeid(Kokkos::Compat::KokkosOpenMPWrapperNode).name())
273 # ifdef HAVE_MUELU_CUDA 274 if (
typeid(
Node).name() ==
typeid(Kokkos::Compat::KokkosCudaWrapperNode).name())
286 if (paramList.isParameter(
"cycle type")) {
287 std::map<std::string, CycleType> cycleMap;
291 auto cycleType = paramList.get<std::string>(
"cycle type");
293 "Invalid cycle type: \"" << cycleType <<
"\"");
294 Cycle_ = cycleMap[cycleType];
297 if (paramList.isParameter(
"W cycle start level")) {
298 WCycleStartLevel_ = paramList.get<
int>(
"W cycle start level");
301 if (paramList.isParameter(
"coarse grid correction scaling factor"))
302 scalingFactor_ = paramList.get<
double>(
"coarse grid correction scaling factor");
304 this->maxCoarseSize_ = paramList.get<
int> (
"coarse: max size", MasterList::getDefault<int>(
"coarse: max size"));
305 this->numDesiredLevel_ = paramList.get<
int> (
"max levels", MasterList::getDefault<int>(
"max levels"));
306 blockSize_ = paramList.get<
int> (
"number of equations", MasterList::getDefault<int>(
"number of equations"));
311 if (paramList.isSublist(
"export data")) {
312 ParameterList printList = paramList.sublist(
"export data");
314 if (printList.isParameter(
"A"))
315 this->matricesToPrint_ = Teuchos::getArrayFromStringParameter<int>(printList,
"A");
316 if (printList.isParameter(
"P"))
317 this->prolongatorsToPrint_ = Teuchos::getArrayFromStringParameter<int>(printList,
"P");
318 if (printList.isParameter(
"R"))
319 this->restrictorsToPrint_ = Teuchos::getArrayFromStringParameter<int>(printList,
"R");
320 if (printList.isParameter(
"Nullspace"))
321 this->nullspaceToPrint_ = Teuchos::getArrayFromStringParameter<int>(printList,
"Nullspace");
322 if (printList.isParameter(
"Coordinates"))
323 this->coordinatesToPrint_ = Teuchos::getArrayFromStringParameter<int>(printList,
"Coordinates");
324 if (printList.isParameter(
"pcoarsen: element to node map"))
325 this->elementToNodeMapsToPrint_ = Teuchos::getArrayFromStringParameter<int>(printList,
"pcoarsen: element to node map");
337 if (outputFilename !=
"")
347 useCoordinates_ =
false;
348 if (
MUELU_TEST_PARAM_2LIST(paramList, paramList,
"aggregation: drop scheme", std::string,
"distance laplacian") ||
351 useCoordinates_ =
true;
352 }
else if(paramList.isSublist(
"smoother: params")) {
353 const auto smooParamList = paramList.sublist(
"smoother: params");
354 if(smooParamList.isParameter(
"partitioner: type") &&
355 (smooParamList.get<std::string>(
"partitioner: type") ==
"line")) {
356 useCoordinates_ =
true;
359 for (
int levelID = 0; levelID < this->numDesiredLevel_; levelID++) {
360 std::string levelStr =
"level " +
toString(levelID);
362 if (paramList.isSublist(levelStr)) {
363 const ParameterList& levelList = paramList.sublist(levelStr);
365 if (
MUELU_TEST_PARAM_2LIST(levelList, paramList,
"aggregation: drop scheme", std::string,
"distance laplacian") ||
368 useCoordinates_ =
true;
376 if (!paramList.isSublist(
"repartition: params")) {
377 useCoordinates_ =
true;
379 const ParameterList& repParams = paramList.sublist(
"repartition: params");
380 if (repParams.isType<std::string>(
"algorithm")) {
381 const std::string algo = repParams.get<std::string>(
"algorithm");
382 if (algo ==
"multijagged" || algo ==
"rcb") {
383 useCoordinates_ =
true;
386 useCoordinates_ =
true;
390 for (
int levelID = 0; levelID < this->numDesiredLevel_; levelID++) {
391 std::string levelStr =
"level " +
toString(levelID);
393 if (paramList.isSublist(levelStr)) {
394 const ParameterList& levelList = paramList.sublist(levelStr);
397 if (!levelList.isSublist(
"repartition: params")) {
398 useCoordinates_ =
true;
401 const ParameterList& repParams = levelList.sublist(
"repartition: params");
402 if (repParams.isType<std::string>(
"algorithm")) {
403 const std::string algo = repParams.get<std::string>(
"algorithm");
404 if (algo ==
"multijagged" || algo ==
"rcb"){
405 useCoordinates_ =
true;
409 useCoordinates_ =
true;
418 changedPRrebalance_ =
false;
420 changedPRrebalance_ =
MUELU_TEST_AND_SET_VAR(paramList,
"repartition: rebalance P and R",
bool, this->doPRrebalance_);
423 changedImplicitTranspose_ =
MUELU_TEST_AND_SET_VAR(paramList,
"transpose: use implicit",
bool, this->implicitTranspose_);
426 MUELU_TEST_AND_SET_VAR(paramList,
"fuse prolongation and update",
bool, this->fuseProlongationAndUpdate_);
428 if (paramList.isSublist(
"matvec params"))
429 this->matvecParams_ = Teuchos::parameterList(paramList.sublist(
"matvec params"));
434 defaultManager->SetVerbLevel(this->verbosity_);
435 defaultManager->SetKokkosRefactor(useKokkos_);
438 std::vector<keep_pair> keeps0;
439 UpdateFactoryManager(paramList, ParameterList(), *defaultManager, 0, keeps0);
442 for (
int levelID = 0; levelID < this->numDesiredLevel_; levelID++) {
448 RCP<FactoryManager> levelManager = rcp(
new FactoryManager(*defaultManager));
449 levelManager->SetVerbLevel(defaultManager->GetVerbLevel());
451 std::vector<keep_pair> keeps;
452 if (paramList.isSublist(
"level " +
toString(levelID))) {
454 ParameterList& levelList = paramList.sublist(
"level " +
toString(levelID),
true);
455 UpdateFactoryManager(levelList, paramList, *levelManager, levelID, keeps);
458 ParameterList levelList;
459 UpdateFactoryManager(levelList, paramList, *levelManager, levelID, keeps);
462 this->keep_[levelID] = keeps;
463 this->AddFactoryManager(levelID, 1, levelManager);
470 this->GetOStream(static_cast<MsgType>(
Runtime1), 0) << paramList << std::endl;
474 ParameterList unusedParamList;
477 for (ParameterList::ConstIterator it = paramList.begin(); it != paramList.end(); it++) {
478 const ParameterEntry& entry = paramList.entry(it);
480 if (!entry.isList() && !entry.isUsed())
481 unusedParamList.setEntry(paramList.name(it), entry);
485 for (
int levelID = 0; levelID < this->numDesiredLevel_; levelID++) {
486 std::string levelStr =
"level " +
toString(levelID);
488 if (paramList.isSublist(levelStr)) {
489 const ParameterList& levelList = paramList.sublist(levelStr);
491 for (ParameterList::ConstIterator itr = levelList.begin(); itr != levelList.end(); ++itr) {
492 const ParameterEntry& entry = levelList.entry(itr);
494 if (!entry.isList() && !entry.isUsed())
495 unusedParamList.sublist(levelStr).setEntry(levelList.name(itr), entry);
500 if (unusedParamList.numParams() > 0) {
501 std::ostringstream unusedParamsStream;
503 unusedParamList.print(unusedParamsStream, indent);
505 this->GetOStream(
Warnings1) <<
"The following parameters were not used:\n" << unusedParamsStream.str() << std::endl;
517 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
520 int levelID, std::vector<keep_pair>& keeps)
const 525 using strings = std::unordered_set<std::string>;
528 if (paramList.numParams() == 0 && defaultList.numParams() > 0)
529 paramList = ParameterList(defaultList);
532 TEUCHOS_TEST_FOR_EXCEPTION(strings({
"none",
"tP",
"RP",
"emin",
"RAP",
"full",
"S"}).count(reuseType) == 0,
535 MUELU_SET_VAR_2LIST(paramList, defaultList,
"multigrid algorithm", std::string, multigridAlgo);
536 TEUCHOS_TEST_FOR_EXCEPTION(strings({
"unsmoothed",
"sa",
"pg",
"emin",
"matlab",
"pcoarsen"}).count(multigridAlgo) == 0,
537 Exceptions::RuntimeError,
"Unknown \"multigrid algorithm\" value: \"" << multigridAlgo <<
"\". Please consult User's Guide.");
538 #ifndef HAVE_MUELU_MATLAB 540 "Cannot use matlab for multigrid algorithm - MueLu was not configured with MATLAB support.");
542 #ifndef HAVE_MUELU_INTREPID2 544 "Cannot use IntrepidPCoarsen prolongator factory - MueLu was not configured with Intrepid support.");
549 if (reuseType ==
"none" || reuseType ==
"S" || reuseType ==
"RP" || reuseType ==
"RAP") {
552 }
else if (reuseType ==
"tP" && (multigridAlgo !=
"sa" && multigridAlgo !=
"unsmoothed")) {
554 this->GetOStream(
Warnings0) <<
"Ignoring \"tP\" reuse option as it is only compatible with \"sa\", " 555 "or \"unsmoothed\" multigrid algorithms" << std::endl;
557 }
else if (reuseType ==
"emin" && multigridAlgo !=
"emin") {
559 this->GetOStream(
Warnings0) <<
"Ignoring \"emin\" reuse option it is only compatible with " 560 "\"emin\" multigrid algorithm" << std::endl;
565 bool have_userP =
false;
566 if (paramList.isParameter(
"P") && !paramList.get<RCP<Matrix> >(
"P").is_null())
570 UpdateFactoryManager_Smoothers(paramList, defaultList, manager, levelID, keeps);
573 UpdateFactoryManager_CoarseSolvers(paramList, defaultList, manager, levelID, keeps);
576 UpdateFactoryManager_Aggregation_TentativeP(paramList, defaultList, manager, levelID, keeps);
579 RCP<Factory> nullSpaceFactory;
580 UpdateFactoryManager_Nullspace(paramList, defaultList, manager, levelID, keeps, nullSpaceFactory);
591 }
else if (multigridAlgo ==
"unsmoothed") {
595 }
else if (multigridAlgo ==
"sa") {
597 UpdateFactoryManager_SA(paramList, defaultList, manager, levelID, keeps);
599 }
else if (multigridAlgo ==
"emin") {
601 UpdateFactoryManager_Emin(paramList, defaultList, manager, levelID, keeps);
603 }
else if (multigridAlgo ==
"pg") {
605 UpdateFactoryManager_PG(paramList, defaultList, manager, levelID, keeps);
607 }
else if (multigridAlgo ==
"matlab") {
609 UpdateFactoryManager_Matlab(paramList, defaultList, manager, levelID, keeps);
611 }
else if (multigridAlgo ==
"pcoarsen") {
613 UpdateFactoryManager_PCoarsen(paramList, defaultList, manager, levelID, keeps);
617 UpdateFactoryManager_SemiCoarsen(paramList, defaultList, manager, levelID, keeps);
620 UpdateFactoryManager_Restriction(paramList, defaultList, manager, levelID, keeps);
623 UpdateFactoryManager_RAP(paramList, defaultList, manager, levelID, keeps);
626 UpdateFactoryManager_Coordinates(paramList, defaultList, manager, levelID, keeps);
629 if ((reuseType ==
"RP" || reuseType ==
"RAP" || reuseType ==
"full") && levelID)
632 if (reuseType ==
"RP" && levelID) {
634 if (!this->implicitTranspose_)
637 if ((reuseType ==
"tP" || reuseType ==
"RP" || reuseType ==
"emin") && useCoordinates_ && levelID)
641 UpdateFactoryManager_Repartition(paramList, defaultList, manager, levelID, keeps, nullSpaceFactory);
644 if ((reuseType ==
"RAP" || reuseType ==
"full") && levelID) {
646 if (!this->implicitTranspose_)
655 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
658 FactoryManager& manager,
int levelID, std::vector<keep_pair>& keeps)
const 660 MUELU_SET_VAR_2LIST(paramList, defaultList,
"multigrid algorithm", std::string, multigridAlgo);
665 bool isCustomSmoother =
666 paramList.isParameter(
"smoother: pre or post") ||
667 paramList.isParameter(
"smoother: type") || paramList.isParameter(
"smoother: pre type") || paramList.isParameter(
"smoother: post type") ||
668 paramList.isSublist (
"smoother: params") || paramList.isSublist (
"smoother: pre params") || paramList.isSublist (
"smoother: post params") ||
669 paramList.isParameter(
"smoother: sweeps") || paramList.isParameter(
"smoother: pre sweeps") || paramList.isParameter(
"smoother: post sweeps") ||
670 paramList.isParameter(
"smoother: overlap") || paramList.isParameter(
"smoother: pre overlap") || paramList.isParameter(
"smoother: post overlap");
674 manager.
SetFactory(
"Smoother", Teuchos::null);
676 }
else if (isCustomSmoother) {
680 #define TEST_MUTUALLY_EXCLUSIVE(arg1,arg2) \ 681 TEUCHOS_TEST_FOR_EXCEPTION(paramList.isParameter(#arg1) && paramList.isParameter(#arg2), \ 682 Exceptions::InvalidArgument, "You cannot specify both \""#arg1"\" and \""#arg2"\""); 683 #define TEST_MUTUALLY_EXCLUSIVE_S(arg1,arg2) \ 684 TEUCHOS_TEST_FOR_EXCEPTION(paramList.isSublist(#arg1) && paramList.isSublist(#arg2), \ 685 Exceptions::InvalidArgument, "You cannot specify both \""#arg1"\" and \""#arg2"\""); 695 TEUCHOS_TEST_FOR_EXCEPTION(
PreOrPost ==
"both" && (paramList.isParameter(
"smoother: pre type") != paramList.isParameter(
"smoother: post type")),
700 ParameterList defaultSmootherParams;
701 defaultSmootherParams.set(
"relaxation: type",
"Symmetric Gauss-Seidel");
702 defaultSmootherParams.set(
"relaxation: sweeps", Teuchos::OrdinalTraits<LO>::one());
703 defaultSmootherParams.set(
"relaxation: damping factor", Teuchos::ScalarTraits<Scalar>::one());
705 RCP<SmootherFactory> preSmoother = Teuchos::null, postSmoother = Teuchos::null;
706 std::string preSmootherType, postSmootherType;
707 ParameterList preSmootherParams, postSmootherParams;
709 if (paramList.isParameter(
"smoother: overlap"))
710 overlap = paramList.get<
int>(
"smoother: overlap");
713 if (paramList.isParameter(
"smoother: pre type")) {
714 preSmootherType = paramList.get<std::string>(
"smoother: pre type");
716 MUELU_SET_VAR_2LIST(paramList, defaultList,
"smoother: type", std::string, preSmootherTypeTmp);
717 preSmootherType = preSmootherTypeTmp;
719 if (paramList.isParameter(
"smoother: pre overlap"))
720 overlap = paramList.get<
int>(
"smoother: pre overlap");
722 if (paramList.isSublist(
"smoother: pre params"))
723 preSmootherParams = paramList.sublist(
"smoother: pre params");
724 else if (paramList.isSublist(
"smoother: params"))
725 preSmootherParams = paramList.sublist(
"smoother: params");
726 else if (defaultList.isSublist(
"smoother: params"))
727 preSmootherParams = defaultList.sublist(
"smoother: params");
728 else if (preSmootherType ==
"RELAXATION")
729 preSmootherParams = defaultSmootherParams;
731 #ifdef HAVE_MUELU_INTREPID2 733 if (multigridAlgo ==
"pcoarsen" && preSmootherType ==
"TOPOLOGICAL" &&
734 defaultList.isParameter(
"pcoarsen: schedule") && defaultList.isParameter(
"pcoarsen: element")) {
737 auto pcoarsen_schedule = Teuchos::getArrayFromStringParameter<int>(defaultList,
"pcoarsen: schedule");
738 auto pcoarsen_element = defaultList.get<std::string>(
"pcoarsen: element");
740 if (levelID < (
int)pcoarsen_schedule.size()) {
742 auto lo = pcoarsen_element + std::to_string(pcoarsen_schedule[levelID]);
743 preSmootherParams.set(
"pcoarsen: hi basis", lo);
748 #ifdef HAVE_MUELU_MATLAB 749 if (preSmootherType ==
"matlab")
757 if (paramList.isParameter(
"smoother: post type"))
758 postSmootherType = paramList.get<std::string>(
"smoother: post type");
760 MUELU_SET_VAR_2LIST(paramList, defaultList,
"smoother: type", std::string, postSmootherTypeTmp);
761 postSmootherType = postSmootherTypeTmp;
764 if (paramList.isSublist(
"smoother: post params"))
765 postSmootherParams = paramList.sublist(
"smoother: post params");
766 else if (paramList.isSublist(
"smoother: params"))
767 postSmootherParams = paramList.sublist(
"smoother: params");
768 else if (defaultList.isSublist(
"smoother: params"))
769 postSmootherParams = defaultList.sublist(
"smoother: params");
770 else if (postSmootherType ==
"RELAXATION")
771 postSmootherParams = defaultSmootherParams;
772 if (paramList.isParameter(
"smoother: post overlap"))
773 overlap = paramList.get<
int>(
"smoother: post overlap");
775 if (postSmootherType == preSmootherType &&
areSame(preSmootherParams, postSmootherParams))
776 postSmoother = preSmoother;
778 #ifdef HAVE_MUELU_INTREPID2 780 if (multigridAlgo ==
"pcoarsen" && preSmootherType ==
"TOPOLOGICAL" &&
781 defaultList.isParameter(
"pcoarsen: schedule") && defaultList.isParameter(
"pcoarsen: element")) {
784 auto pcoarsen_schedule = Teuchos::getArrayFromStringParameter<int>(defaultList,
"pcoarsen: schedule");
785 auto pcoarsen_element = defaultList.get<std::string>(
"pcoarsen: element");
787 if (levelID < (
int)pcoarsen_schedule.size()) {
789 auto lo = pcoarsen_element + std::to_string(pcoarsen_schedule[levelID]);
790 postSmootherParams.set(
"pcoarsen: hi basis", lo);
795 #ifdef HAVE_MUELU_MATLAB 796 if (postSmootherType ==
"matlab")
804 if (preSmoother == postSmoother)
807 manager.
SetFactory(
"PreSmoother", preSmoother);
808 manager.
SetFactory(
"PostSmoother", postSmoother);
815 bool reuseSmoothers = (reuseType ==
"S" || reuseType !=
"none");
816 if (reuseSmoothers) {
817 auto preSmootherFactory = rcp_const_cast<
Factory>(rcp_dynamic_cast<
const Factory>(manager.
GetFactory(
"PreSmoother")));
819 if (preSmootherFactory != Teuchos::null) {
820 ParameterList postSmootherFactoryParams;
821 postSmootherFactoryParams.set(
"keep smoother data",
true);
822 preSmootherFactory->SetParameterList(postSmootherFactoryParams);
824 keeps.push_back(
keep_pair(
"PreSmoother data", preSmootherFactory.get()));
827 auto postSmootherFactory = rcp_const_cast<
Factory>(rcp_dynamic_cast<
const Factory>(manager.
GetFactory(
"PostSmoother")));
828 if (postSmootherFactory != Teuchos::null) {
829 ParameterList postSmootherFactoryParams;
830 postSmootherFactoryParams.set(
"keep smoother data",
true);
831 postSmootherFactory->SetParameterList(postSmootherFactoryParams);
833 keeps.push_back(
keep_pair(
"PostSmoother data", postSmootherFactory.get()));
837 if (coarseFactory != Teuchos::null) {
838 ParameterList coarseFactoryParams;
839 coarseFactoryParams.set(
"keep smoother data",
true);
840 coarseFactory->SetParameterList(coarseFactoryParams);
842 keeps.push_back(
keep_pair(
"PreSmoother data", coarseFactory.get()));
846 if ((reuseType ==
"RAP" && levelID) || (reuseType ==
"full")) {
865 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
871 bool isCustomCoarseSolver =
872 paramList.isParameter(
"coarse: type") ||
873 paramList.isParameter(
"coarse: params");
875 this->GetOStream(
Warnings0) <<
"No coarse grid solver" << std::endl;
876 manager.
SetFactory(
"CoarseSolver", Teuchos::null);
878 }
else if (isCustomCoarseSolver) {
885 if (paramList.isParameter(
"coarse: overlap"))
886 overlap = paramList.get<
int>(
"coarse: overlap");
888 ParameterList coarseParams;
889 if (paramList.isSublist(
"coarse: params"))
890 coarseParams = paramList.sublist(
"coarse: params");
891 else if (defaultList.isSublist(
"coarse: params"))
892 coarseParams = defaultList.sublist(
"coarse: params");
894 using strings = std::unordered_set<std::string>;
896 RCP<SmootherPrototype> coarseSmoother;
900 if (strings({
"RELAXATION",
"CHEBYSHEV",
"ILUT",
"ILU",
"RILUK",
"SCHWARZ",
"Amesos",
901 "BLOCK RELAXATION",
"BLOCK_RELAXATION",
"BLOCKRELAXATION" ,
902 "SPARSE BLOCK RELAXATION",
"SPARSE_BLOCK_RELAXATION",
"SPARSEBLOCKRELAXATION",
903 "LINESMOOTHING_BANDEDRELAXATION",
"LINESMOOTHING_BANDED_RELAXATION",
"LINESMOOTHING_BANDED RELAXATION",
904 "LINESMOOTHING_TRIDIRELAXATION",
"LINESMOOTHING_TRIDI_RELAXATION",
"LINESMOOTHING_TRIDI RELAXATION",
905 "LINESMOOTHING_TRIDIAGONALRELAXATION",
"LINESMOOTHING_TRIDIAGONAL_RELAXATION",
"LINESMOOTHING_TRIDIAGONAL RELAXATION",
906 "TOPOLOGICAL",
"FAST_ILU",
"FAST_IC",
"FAST_ILDL"}).count(coarseType)) {
907 coarseSmoother = rcp(
new TrilinosSmoother(coarseType, coarseParams, overlap));
909 #ifdef HAVE_MUELU_MATLAB 910 if (coarseType ==
"matlab")
914 coarseSmoother = rcp(
new DirectSolver(coarseType, coarseParams));
924 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
927 FactoryManager& manager,
int levelID, std::vector<keep_pair>& keeps)
const 929 using strings = std::unordered_set<std::string>;
934 RCP<Factory> dropFactory;
937 #ifdef HAVE_MUELU_MATLAB 939 ParameterList socParams = paramList.sublist(
"strength-of-connection: params");
940 dropFactory->SetParameterList(socParams);
942 throw std::runtime_error(
"Cannot use MATLAB evolutionary strength-of-connection - MueLu was not configured with MATLAB support.");
944 }
else if (
MUELU_TEST_PARAM_2LIST(paramList, paramList,
"aggregation: drop scheme", std::string,
"unsupported vector smoothing")) {
946 ParameterList dropParams;
951 dropFactory->SetParameterList(dropParams);
955 ParameterList dropParams;
956 dropParams.set(
"lightweight wrap",
true);
967 dropFactory->SetParameterList(dropParams);
973 TEUCHOS_TEST_FOR_EXCEPTION(!strings({
"uncoupled",
"coupled",
"brick",
"matlab",
"notay"}).count(aggType),
975 #ifndef HAVE_MUELU_MATLAB 976 if (aggType ==
"matlab")
977 throw std::runtime_error(
"Cannot use MATLAB aggregation - MueLu was not configured with MATLAB support.");
979 RCP<Factory> aggFactory;
980 if (aggType ==
"uncoupled") {
982 ParameterList aggParams;
1003 aggFactory->SetParameterList(aggParams);
1005 aggFactory->SetFactory(
"DofsPerNode", manager.
GetFactory(
"Graph"));
1006 aggFactory->SetFactory(
"Graph", manager.
GetFactory(
"Graph"));
1008 }
else if (aggType ==
"coupled") {
1010 aggFactory->SetFactory(
"Graph", manager.
GetFactory(
"Graph"));
1012 }
else if (aggType ==
"brick") {
1014 ParameterList aggParams;
1021 aggFactory->SetParameterList(aggParams);
1027 aggFactory->SetFactory(
"Coordinates", this->GetFactoryManager(levelID-1)->GetFactory(
"Coordinates"));
1030 else if (aggType ==
"notay") {
1032 ParameterList aggParams;
1038 aggFactory->SetParameterList(aggParams);
1039 aggFactory->SetFactory(
"DofsPerNode", manager.
GetFactory(
"Graph"));
1040 aggFactory->SetFactory(
"Graph", manager.
GetFactory(
"Graph"));
1042 #ifdef HAVE_MUELU_MATLAB 1043 else if(aggType ==
"matlab") {
1044 ParameterList aggParams = paramList.sublist(
"aggregation: params");
1046 aggFactory->SetParameterList(aggParams);
1051 manager.
SetFactory(
"Aggregates", aggFactory);
1055 coarseMap->SetFactory(
"Aggregates", manager.
GetFactory(
"Aggregates"));
1059 if (
MUELU_TEST_PARAM_2LIST(paramList, defaultList,
"aggregation: compute aggregate qualities",
bool,
true)) {
1061 ParameterList aggQualityParams;
1070 aggQualityFact->SetParameterList(aggQualityParams);
1071 manager.
SetFactory(
"AggregateQualities", aggQualityFact);
1073 assert(aggType ==
"uncoupled");
1074 aggFactory->SetFactory(
"AggregateQualities", aggQualityFact);
1080 ParameterList ptentParams;
1081 if (paramList.isSublist(
"matrixmatrix: kernel params"))
1082 ptentParams.sublist(
"matrixmatrix: kernel params",
false) = paramList.sublist(
"matrixmatrix: kernel params");
1083 if (defaultList.isSublist(
"matrixmatrix: kernel params"))
1084 ptentParams.sublist(
"matrixmatrix: kernel params",
false) = defaultList.sublist(
"matrixmatrix: kernel params");
1087 Ptent->SetParameterList(ptentParams);
1088 Ptent->SetFactory(
"Aggregates", manager.
GetFactory(
"Aggregates"));
1089 Ptent->SetFactory(
"CoarseMap", manager.
GetFactory(
"CoarseMap"));
1092 if (reuseType ==
"tP" && levelID) {
1093 keeps.push_back(
keep_pair(
"Nullspace", Ptent.get()));
1094 keeps.push_back(
keep_pair(
"P", Ptent.get()));
1101 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1104 int , std::vector<keep_pair>& keeps)
const 1106 if (paramList.isParameter(
"A") && !paramList.get<RCP<Matrix> >(
"A").is_null()) {
1112 ParameterList RAPparams;
1114 RCP<RAPFactory> RAP;
1115 RCP<RAPShiftFactory> RAPs;
1118 std::string alg = paramList.get(
"rap: algorithm",
"galerkin");
1119 if (alg ==
"shift" || alg ==
"non-galerkin") {
1133 if (paramList.isSublist(
"matrixmatrix: kernel params"))
1134 RAPparams.sublist(
"matrixmatrix: kernel params",
false) = paramList.sublist(
"matrixmatrix: kernel params");
1135 if (defaultList.isSublist(
"matrixmatrix: kernel params"))
1136 RAPparams.sublist(
"matrixmatrix: kernel params",
false) = defaultList.sublist(
"matrixmatrix: kernel params");
1141 if (!paramList.isParameter(
"rap: triple product") &&
1142 paramList.isType<std::string>(
"multigrid algorithm") &&
1143 paramList.get<std::string>(
"multigrid algorithm") ==
"unsmoothed")
1144 paramList.set(
"rap: triple product",
true);
1149 if (paramList.isParameter(
"aggregation: allow empty prolongator columns")) {
1150 RAPparams.set(
"CheckMainDiagonal", paramList.get<
bool>(
"aggregation: allow empty prolongator columns"));
1151 RAPparams.set(
"RepairMainDiagonal", paramList.get<
bool>(
"aggregation: allow empty prolongator columns"));
1153 else if (defaultList.isParameter(
"aggregation: allow empty prolongator columns")) {
1154 RAPparams.set(
"CheckMainDiagonal", defaultList.get<
bool>(
"aggregation: allow empty prolongator columns"));
1155 RAPparams.set(
"RepairMainDiagonal", defaultList.get<
bool>(
"aggregation: allow empty prolongator columns"));
1158 }
catch (Teuchos::Exceptions::InvalidParameterType&) {
1159 TEUCHOS_TEST_FOR_EXCEPTION_PURE_MSG(
true, Teuchos::Exceptions::InvalidParameterType,
1160 "Error: parameter \"aggregation: allow empty prolongator columns\" must be of type " << Teuchos::TypeNameTraits<bool>::name());
1163 if (!RAP.is_null()) {
1164 RAP->SetParameterList(RAPparams);
1165 RAP->SetFactory(
"P", manager.
GetFactory(
"P"));
1167 RAPs->SetParameterList(RAPparams);
1168 RAPs->SetFactory(
"P", manager.
GetFactory(
"P"));
1171 if (!this->implicitTranspose_) {
1173 RAP->SetFactory(
"R", manager.
GetFactory(
"R"));
1175 RAPs->SetFactory(
"R", manager.
GetFactory(
"R"));
1178 if (
MUELU_TEST_PARAM_2LIST(paramList, defaultList,
"aggregation: export visualization data",
bool,
true)) {
1180 ParameterList aggExportParams;
1188 aggExport->SetParameterList(aggExportParams);
1189 aggExport->SetFactory(
"DofsPerNode", manager.
GetFactory(
"DofsPerNode"));
1192 RAP->AddTransferFactory(aggExport);
1194 RAPs->AddTransferFactory(aggExport);
1202 MUELU_SET_VAR_2LIST(paramList, defaultList,
"sa: use filtered matrix",
bool, useFiltering);
1203 bool filteringChangesMatrix = useFiltering && !
MUELU_TEST_PARAM_2LIST(paramList, defaultList,
"aggregation: drop tol",
double, 0);
1205 if (reuseType ==
"RP" || (reuseType ==
"tP" && !filteringChangesMatrix)) {
1206 if (!RAP.is_null()) {
1207 keeps.push_back(
keep_pair(
"AP reuse data", RAP.get()));
1208 keeps.push_back(
keep_pair(
"RAP reuse data", RAP.get()));
1211 keeps.push_back(
keep_pair(
"AP reuse data", RAPs.get()));
1212 keeps.push_back(
keep_pair(
"RAP reuse data", RAPs.get()));
1220 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1225 bool have_userCO =
false;
1226 if (paramList.isParameter(
"Coordinates") && !paramList.get<RCP<MultiVector> >(
"Coordinates").is_null())
1229 if (useCoordinates_) {
1235 coords->SetFactory(
"Aggregates", manager.
GetFactory(
"Aggregates"));
1236 coords->SetFactory(
"CoarseMap", manager.
GetFactory(
"CoarseMap"));
1240 if (!RAP.is_null()) {
1241 RAP->AddTransferFactory(manager.
GetFactory(
"Coordinates"));
1244 RAPs->AddTransferFactory(manager.
GetFactory(
"Coordinates"));
1253 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1256 int levelID, std::vector<keep_pair>& )
const 1258 MUELU_SET_VAR_2LIST(paramList, defaultList,
"multigrid algorithm", std::string, multigridAlgo);
1259 bool have_userR =
false;
1260 if (paramList.isParameter(
"R") && !paramList.get<RCP<Matrix> >(
"R").is_null())
1265 if (!this->implicitTranspose_) {
1268 if (isSymmetric ==
false && (multigridAlgo ==
"unsmoothed" || multigridAlgo ==
"emin")) {
1270 "Switching \"problem: symmetric\" parameter to symmetric as multigrid algorithm. " <<
1271 multigridAlgo <<
" is primarily supposed to be used for symmetric problems.\n\n" <<
1272 "Please note: if you are using \"unsmoothed\" transfer operators the \"problem: symmetric\" parameter " <<
1273 "has no real mathematical meaning, i.e. you can use it for non-symmetric\n" <<
1274 "problems, too. With \"problem: symmetric\"=\"symmetric\" you can use implicit transpose for building " <<
1275 "the restriction operators which may drastically reduce the amount of consumed memory." << std::endl;
1279 "Petrov-Galerkin smoothed transfer operators are only allowed for non-symmetric problems: Set \"problem: symmetric\" to false!\n" \
1280 "While PG smoothed transfer operators generally would also work for symmetric problems this is an unusual use case. " \
1281 "You can use the factory-based xml interface though if you need PG-AMG for symmetric problems.");
1298 if (paramList.isParameter(
"restriction: scale nullspace") && paramList.get<
bool>(
"restriction: scale nullspace")) {
1300 Teuchos::ParameterList tentPlist;
1301 tentPlist.set(
"Nullspace name",
"Scaled Nullspace");
1302 tentPFactory->SetParameterList(tentPlist);
1303 tentPFactory->SetFactory(
"Aggregates",manager.
GetFactory(
"Aggregates"));
1304 tentPFactory->SetFactory(
"CoarseMap",manager.
GetFactory(
"CoarseMap"));
1307 R->SetFactory(
"P",tentPFactory);
1316 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1319 int levelID, std::vector<keep_pair>& keeps, RCP<Factory> & nullSpaceFactory)
const 1324 MUELU_SET_VAR_2LIST(paramList, defaultList,
"repartition: node repartition level",
int,nodeRepartitionLevel);
1364 "Reuse types \"tP\" and \"PR\" require \"repartition: rebalance P and R\" set to \"false\"");
1369 MUELU_SET_VAR_2LIST(paramList, defaultList,
"repartition: partitioner", std::string, partName);
1371 "Invalid partitioner name: \"" << partName <<
"\". Valid options: \"zoltan\", \"zoltan2\"");
1373 #ifndef HAVE_MUELU_ZOLTAN 1374 bool switched =
false;
1375 if (partName ==
"zoltan") {
1376 this->GetOStream(
Warnings0) <<
"Zoltan interface is not available, trying to switch to Zoltan2" << std::endl;
1377 partName =
"zoltan2";
1381 # ifndef HAVE_MUELU_ZOLTAN2 1382 bool switched =
false;
1385 #ifndef HAVE_MUELU_ZOLTAN2 1386 if (partName ==
"zoltan2" && !switched) {
1387 this->GetOStream(
Warnings0) <<
"Zoltan2 interface is not available, trying to switch to Zoltan" << std::endl;
1388 partName =
"zoltan";
1394 ParameterList repartheurParams;
1400 repartheurFactory->SetParameterList(repartheurParams);
1401 repartheurFactory->SetFactory(
"A", manager.
GetFactory(
"A"));
1402 manager.
SetFactory(
"number of partitions", repartheurFactory);
1403 manager.
SetFactory(
"repartition: heuristic target rows per process", repartheurFactory);
1406 RCP<Factory> partitioner;
1407 if (levelID == nodeRepartitionLevel) {
1411 ParameterList partParams;
1413 partitioner->SetParameterList(partParams);
1414 partitioner->SetFactory(
"Node Comm", manager.
GetFactory(
"Node Comm"));
1419 else if (partName ==
"zoltan") {
1420 #ifdef HAVE_MUELU_ZOLTAN 1426 }
else if (partName ==
"zoltan2") {
1427 #ifdef HAVE_MUELU_ZOLTAN2 1429 ParameterList partParams;
1430 RCP<const ParameterList> partpartParams = rcp(
new ParameterList(paramList.sublist(
"repartition: params",
false)));
1431 partParams.set(
"ParameterList", partpartParams);
1432 partitioner->SetParameterList(partParams);
1433 partitioner->SetFactory(
"repartition: heuristic target rows per process",
1434 manager.
GetFactory(
"repartition: heuristic target rows per process"));
1440 partitioner->SetFactory(
"A", manager.
GetFactory(
"A"));
1441 partitioner->SetFactory(
"number of partitions", manager.
GetFactory(
"number of partitions"));
1442 if (useCoordinates_)
1443 partitioner->SetFactory(
"Coordinates", manager.
GetFactory(
"Coordinates"));
1444 manager.
SetFactory(
"Partition", partitioner);
1448 ParameterList repartParams;
1452 repartFactory->SetParameterList(repartParams);
1453 repartFactory->SetFactory(
"A", manager.
GetFactory(
"A"));
1454 repartFactory->SetFactory(
"number of partitions", manager.
GetFactory(
"number of partitions"));
1455 repartFactory->SetFactory(
"Partition", manager.
GetFactory(
"Partition"));
1456 manager.
SetFactory(
"Importer", repartFactory);
1457 if (reuseType !=
"none" && reuseType !=
"S" && levelID)
1462 ParameterList rebAcParams;
1464 newA->SetParameterList(rebAcParams);
1465 newA->SetFactory(
"A", manager.
GetFactory(
"A"));
1466 newA->SetFactory(
"Importer", manager.
GetFactory(
"Importer"));
1471 ParameterList newPparams;
1472 newPparams.set(
"type",
"Interpolation");
1473 if (changedPRrebalance_)
1474 newPparams.set(
"repartition: rebalance P and R", this->doPRrebalance_);
1476 newP-> SetParameterList(newPparams);
1477 newP-> SetFactory(
"Importer", manager.
GetFactory(
"Importer"));
1478 newP-> SetFactory(
"P", manager.
GetFactory(
"P"));
1479 if (!paramList.isParameter(
"semicoarsen: number of levels"))
1480 newP->SetFactory(
"Nullspace", manager.
GetFactory(
"Ptent"));
1482 newP->SetFactory(
"Nullspace", manager.
GetFactory(
"P"));
1483 if (useCoordinates_)
1484 newP-> SetFactory(
"Coordinates", manager.
GetFactory(
"Coordinates"));
1486 if (useCoordinates_)
1491 ParameterList newRparams;
1492 newRparams.set(
"type",
"Restriction");
1494 if (changedPRrebalance_)
1495 newRparams.set(
"repartition: rebalance P and R", this->doPRrebalance_);
1496 if (changedImplicitTranspose_)
1497 newRparams.set(
"transpose: use implicit", this->implicitTranspose_);
1498 newR-> SetParameterList(newRparams);
1499 newR-> SetFactory(
"Importer", manager.
GetFactory(
"Importer"));
1500 if (!this->implicitTranspose_) {
1501 newR->SetFactory(
"R", manager.
GetFactory(
"R"));
1512 nullSpaceFactory->SetFactory(
"Nullspace", newP);
1522 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1525 int , std::vector<keep_pair>& , RCP<Factory> & nullSpaceFactory)
const 1530 bool have_userNS =
false;
1531 if (paramList.isParameter(
"Nullspace") && !paramList.get<RCP<MultiVector> >(
"Nullspace").is_null())
1535 nullSpace->SetFactory(
"Nullspace", manager.
GetFactory(
"Ptent"));
1538 nullSpaceFactory = nullSpace;
1540 if (paramList.isParameter(
"restriction: scale nullspace") && paramList.get<
bool>(
"restriction: scale nullspace")) {
1542 scaledNSfactory->SetFactory(
"Nullspace",nullSpaceFactory);
1543 manager.
SetFactory(
"Scaled Nullspace",scaledNSfactory);
1551 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1554 int , std::vector<keep_pair>& )
const 1557 RCP<SemiCoarsenPFactory> semicoarsenFactory = Teuchos::null;
1558 if (paramList.isParameter(
"semicoarsen: number of levels") &&
1559 paramList.get<
int>(
"semicoarsen: number of levels") > 0) {
1561 ParameterList togglePParams;
1562 ParameterList semicoarsenPParams;
1563 ParameterList linedetectionParams;
1573 linedetectionFactory->SetParameterList(linedetectionParams);
1574 semicoarsenFactory ->SetParameterList(semicoarsenPParams);
1575 togglePFactory ->SetParameterList(togglePParams);
1577 togglePFactory->AddCoarseNullspaceFactory (semicoarsenFactory);
1578 togglePFactory->AddProlongatorFactory (semicoarsenFactory);
1579 togglePFactory->AddPtentFactory (semicoarsenFactory);
1580 togglePFactory->AddCoarseNullspaceFactory (manager.
GetFactory(
"Ptent"));
1581 togglePFactory->AddProlongatorFactory (manager.
GetFactory(
"P"));
1582 togglePFactory->AddPtentFactory (manager.
GetFactory(
"Ptent"));
1584 manager.
SetFactory(
"CoarseNumZLayers", linedetectionFactory);
1585 manager.
SetFactory(
"LineDetection_Layers", linedetectionFactory);
1586 manager.
SetFactory(
"LineDetection_VertLineIds", linedetectionFactory);
1590 manager.
SetFactory(
"Nullspace", togglePFactory);
1594 if (paramList.isParameter(
"semicoarsen: number of levels")) {
1596 tf->SetFactory(
"Chosen P", manager.
GetFactory(
"P"));
1597 tf->AddCoordTransferFactory(semicoarsenFactory);
1600 coords->SetFactory(
"Aggregates", manager.
GetFactory(
"Aggregates"));
1601 coords->SetFactory(
"CoarseMap", manager.
GetFactory(
"CoarseMap"));
1602 tf->AddCoordTransferFactory(coords);
1611 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1614 int levelID, std::vector<keep_pair>& keeps)
const 1616 #ifdef HAVE_MUELU_INTREPID2 1618 if (defaultList.isParameter(
"pcoarsen: schedule") && defaultList.isParameter(
"pcoarsen: element")) {
1621 auto pcoarsen_schedule = Teuchos::getArrayFromStringParameter<int>(defaultList,
"pcoarsen: schedule");
1622 auto pcoarsen_element = defaultList.get<std::string>(
"pcoarsen: element");
1624 if (levelID >= (
int)pcoarsen_schedule.size()) {
1627 UpdateFactoryManager_SA(paramList, defaultList, manager, levelID, keeps);
1631 ParameterList Pparams;
1633 std::string lo = pcoarsen_element + std::to_string(pcoarsen_schedule[levelID]);
1634 std::string hi = (levelID ? pcoarsen_element + std::to_string(pcoarsen_schedule[levelID-1]) : lo);
1635 Pparams.set(
"pcoarsen: hi basis", hi);
1636 Pparams.set(
"pcoarsen: lo basis", lo);
1637 P->SetParameterList(Pparams);
1646 ParameterList Pparams;
1650 P->SetParameterList(Pparams);
1663 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1668 ParameterList Pparams;
1669 if (paramList.isSublist(
"matrixmatrix: kernel params"))
1670 Pparams.sublist(
"matrixmatrix: kernel params",
false) = paramList.sublist(
"matrixmatrix: kernel params");
1671 if (defaultList.isSublist(
"matrixmatrix: kernel params"))
1672 Pparams.sublist(
"matrixmatrix: kernel params",
false) = defaultList.sublist(
"matrixmatrix: kernel params");
1674 P->SetParameterList(Pparams);
1677 MUELU_SET_VAR_2LIST(paramList, defaultList,
"sa: use filtered matrix",
bool, useFiltering);
1685 ParameterList fParams;
1689 filterFactory->SetParameterList(fParams);
1690 filterFactory->SetFactory(
"Graph", manager.
GetFactory(
"Graph"));
1692 filterFactory->SetFactory(
"Filtering", manager.
GetFactory(
"Graph"));
1694 P->SetFactory(
"A", filterFactory);
1697 P->SetFactory(
"A", manager.
GetFactory(
"Graph"));
1701 P->SetFactory(
"P", manager.
GetFactory(
"Ptent"));
1704 bool filteringChangesMatrix = useFiltering && !
MUELU_TEST_PARAM_2LIST(paramList, defaultList,
"aggregation: drop tol",
double, 0);
1706 if (reuseType ==
"tP" && !filteringChangesMatrix)
1707 keeps.push_back(
keep_pair(
"AP reuse data", P.get()));
1713 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1716 int , std::vector<keep_pair>& )
const 1721 "Invalid pattern name: \"" << patternType <<
"\". Valid options: \"AkPtent\"");
1724 ParameterList patternParams;
1726 patternFactory->SetParameterList(patternParams);
1727 patternFactory->SetFactory(
"P", manager.
GetFactory(
"Ptent"));
1728 manager.
SetFactory(
"Ppattern", patternFactory);
1732 constraintFactory->SetFactory(
"Ppattern", manager.
GetFactory(
"Ppattern"));
1733 constraintFactory->SetFactory(
"CoarseNullspace", manager.
GetFactory(
"Ptent"));
1734 manager.
SetFactory(
"Constraint", constraintFactory);
1739 MUELU_SET_VAR_2LIST(paramList, defaultList,
"emin: use filtered matrix",
bool, useFiltering);
1747 ParameterList fParams;
1751 filterFactory->SetParameterList(fParams);
1752 filterFactory->SetFactory(
"Graph", manager.
GetFactory(
"Graph"));
1754 filterFactory->SetFactory(
"Filtering", manager.
GetFactory(
"Graph"));
1756 P->SetFactory(
"A", filterFactory);
1759 P->SetFactory(
"A", manager.
GetFactory(
"Graph"));
1764 ParameterList Pparams;
1767 if (reuseType ==
"emin") {
1769 Pparams.set(
"Keep P0",
true);
1770 Pparams.set(
"Keep Constraint0",
true);
1772 P->SetParameterList(Pparams);
1773 P->SetFactory(
"P", manager.
GetFactory(
"Ptent"));
1774 P->SetFactory(
"Constraint", manager.
GetFactory(
"Constraint"));
1781 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1784 int , std::vector<keep_pair>& )
const 1787 "Implicit transpose not supported with Petrov-Galerkin smoothed transfer operators: Set \"transpose: use implicit\" to false!\n" \
1788 "Petrov-Galerkin transfer operator smoothing for non-symmetric problems requires a separate handling of the restriction operator which " \
1789 "does not allow the usage of implicit transpose easily.");
1793 P->SetFactory(
"P", manager.
GetFactory(
"Ptent"));
1801 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1804 int , std::vector<keep_pair>& )
const {
1805 #ifdef HAVE_MUELU_MATLAB 1806 ParameterList Pparams = paramList.sublist(
"transfer: params");
1808 P->SetParameterList(Pparams);
1809 P->SetFactory(
"P", manager.
GetFactory(
"Ptent"));
1817 #undef MUELU_SET_VAR_2LIST 1818 #undef MUELU_TEST_AND_SET_VAR 1819 #undef MUELU_TEST_AND_SET_PARAM_2LIST 1820 #undef MUELU_TEST_PARAM_2LIST 1821 #undef MUELU_KOKKOS_FACTORY 1825 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1827 ParameterList paramList = constParamList;
1830 const int maxLevels = 100;
1833 std::vector<ParameterList> paramLists;
1834 for (
int levelID = 0; levelID < maxLevels; levelID++) {
1835 std::string sublistName =
"level " +
toString(levelID);
1836 if (paramList.isSublist(sublistName)) {
1837 paramLists.push_back(paramList.sublist(sublistName));
1839 paramList.remove(sublistName);
1842 paramLists.push_back(paramList);
1844 #ifdef HAVE_MUELU_MATLAB 1846 for (
size_t i = 0; i < paramLists.size(); i++) {
1847 std::vector<std::string> customVars;
1849 for(Teuchos::ParameterList::ConstIterator it = paramLists[i].begin(); it != paramLists[i].end(); it++) {
1850 std::string paramName = paramLists[i].name(it);
1853 customVars.push_back(paramName);
1857 for (
size_t j = 0; j < customVars.size(); j++)
1858 paramLists[i].
remove(customVars[j],
false);
1862 const int maxDepth = 0;
1863 for (
size_t i = 0; i < paramLists.size(); i++) {
1866 paramLists[i].validateParameters(validList, maxDepth);
1868 }
catch (
const Teuchos::Exceptions::InvalidParameterName& e) {
1869 std::string eString = e.what();
1872 size_t nameStart = eString.find_first_of(
'"') + 1;
1873 size_t nameEnd = eString.find_first_of(
'"', nameStart);
1874 std::string name = eString.substr(nameStart, nameEnd - nameStart);
1876 size_t bestScore = 100;
1877 std::string bestName =
"";
1878 for (ParameterList::ConstIterator it = validList.begin(); it != validList.end(); it++) {
1879 const std::string& pName = validList.name(it);
1880 this->GetOStream(
Runtime1) <<
"| " << pName;
1881 size_t score =
LevenshteinDistance(name.c_str(), name.length(), pName.c_str(), pName.length());
1882 this->GetOStream(
Runtime1) <<
" -> " << score << std::endl;
1883 if (score < bestScore) {
1888 if (bestScore < 10 && bestName !=
"") {
1889 TEUCHOS_TEST_FOR_EXCEPTION(
true, Teuchos::Exceptions::InvalidParameterName,
1890 eString <<
"The parameter name \"" + name +
"\" is not valid. Did you mean \"" + bestName <<
"\"?\n");
1893 TEUCHOS_TEST_FOR_EXCEPTION(
true, Teuchos::Exceptions::InvalidParameterName,
1894 eString <<
"The parameter name \"" + name +
"\" is not valid.\n");
1903 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
1908 ParameterList paramList = constParamList;
1915 if (paramList.isSublist(
"Matrix")) {
1916 blockSize_ = paramList.sublist(
"Matrix").get<
int>(
"PDE equations", MasterList::getDefault<int>(
"number of equations"));
1917 dofOffset_ = paramList.sublist(
"Matrix").get<
GlobalOrdinal>(
"DOF offset", 0);
1921 if (factFact_ == Teuchos::null)
1933 if (paramList.isSublist(
"Factories"))
1934 this->BuildFactoryMap(paramList.sublist(
"Factories"), factoryMap, factoryMap, factoryManagers);
1948 if (paramList.isSublist(
"Hierarchy")) {
1949 ParameterList hieraList = paramList.sublist(
"Hierarchy");
1952 if (hieraList.isParameter(
"max levels")) {
1953 this->numDesiredLevel_ = hieraList.get<
int>(
"max levels");
1954 hieraList.remove(
"max levels");
1957 if (hieraList.isParameter(
"coarse: max size")) {
1958 this->maxCoarseSize_ = hieraList.get<
int>(
"coarse: max size");
1959 hieraList.remove(
"coarse: max size");
1962 if (hieraList.isParameter(
"repartition: rebalance P and R")) {
1963 this->doPRrebalance_ = hieraList.get<
bool>(
"repartition: rebalance P and R");
1964 hieraList.remove(
"repartition: rebalance P and R");
1967 if (hieraList.isParameter(
"transpose: use implicit")) {
1968 this->implicitTranspose_ = hieraList.get<
bool>(
"transpose: use implicit");
1969 hieraList.remove(
"transpose: use implicit");
1972 if (hieraList.isParameter(
"fuse prolongation and update")) {
1973 this->fuseProlongationAndUpdate_ = hieraList.get<
bool>(
"fuse prolongation and update");
1974 hieraList.remove(
"fuse prolongation and update");
1977 if (hieraList.isParameter(
"number of vectors")) {
1978 this->numDesiredLevel_ = hieraList.get<
int>(
"number of vectors");
1979 hieraList.remove(
"number of vectors");
1982 if (hieraList.isSublist(
"matvec params"))
1983 this->matvecParams_ = Teuchos::parameterList(hieraList.sublist(
"matvec params"));
1986 if (hieraList.isParameter(
"coarse grid correction scaling factor")) {
1987 this->scalingFactor_ = hieraList.get<
double>(
"coarse grid correction scaling factor");
1988 hieraList.remove(
"coarse grid correction scaling factor");
1992 if (hieraList.isParameter(
"cycle type")) {
1993 std::map<std::string, CycleType> cycleMap;
1997 std::string cycleType = hieraList.get<std::string>(
"cycle type");
1998 TEUCHOS_TEST_FOR_EXCEPTION(cycleMap.count(cycleType) == 0,
Exceptions::RuntimeError,
"Invalid cycle type: \"" << cycleType <<
"\"");
1999 this->Cycle_ = cycleMap[cycleType];
2002 if (hieraList.isParameter(
"W cycle start level")) {
2003 this->WCycleStartLevel_ = hieraList.get<
int>(
"W cycle start level");
2006 if (hieraList.isParameter(
"verbosity")) {
2007 std::string vl = hieraList.get<std::string>(
"verbosity");
2008 hieraList.remove(
"verbosity");
2012 if (hieraList.isParameter(
"output filename"))
2015 if (hieraList.isParameter(
"dependencyOutputLevel"))
2016 this->graphOutputLevel_ = hieraList.get<
int>(
"dependencyOutputLevel");
2019 if (hieraList.isParameter(
"reuse"))
2022 if (hieraList.isSublist(
"DataToWrite")) {
2025 ParameterList foo = hieraList.sublist(
"DataToWrite");
2026 std::string dataName =
"Matrices";
2027 if (foo.isParameter(dataName))
2028 this->matricesToPrint_ = Teuchos::getArrayFromStringParameter<int>(foo, dataName);
2029 dataName =
"Prolongators";
2030 if (foo.isParameter(dataName))
2031 this->prolongatorsToPrint_ = Teuchos::getArrayFromStringParameter<int>(foo, dataName);
2032 dataName =
"Restrictors";
2033 if (foo.isParameter(dataName))
2034 this->restrictorsToPrint_ = Teuchos::getArrayFromStringParameter<int>(foo, dataName);
2038 for (ParameterList::ConstIterator param = hieraList.begin(); param != hieraList.end(); ++param) {
2039 const std::string & paramName = hieraList.name(param);
2041 if (paramName !=
"DataToWrite" && hieraList.isSublist(paramName)) {
2042 ParameterList levelList = hieraList.sublist(paramName);
2044 int startLevel = 0;
if(levelList.isParameter(
"startLevel")) { startLevel = levelList.get<
int>(
"startLevel"); levelList.remove(
"startLevel"); }
2045 int numDesiredLevel = 1;
if(levelList.isParameter(
"numDesiredLevel")) { numDesiredLevel = levelList.get<
int>(
"numDesiredLevel"); levelList.remove(
"numDesiredLevel"); }
2058 BuildFactoryMap(levelList, factoryMap, levelFactoryMap, factoryManagers);
2060 RCP<FactoryManager> m = rcp(
new FactoryManager(levelFactoryMap));
2061 if (hieraList.isParameter(
"use kokkos refactor"))
2062 m->SetKokkosRefactor(hieraList.get<
bool>(
"use kokkos refactor"));
2064 if (startLevel >= 0)
2065 this->AddFactoryManager(startLevel, numDesiredLevel, m);
2067 TEUCHOS_TEST_FOR_EXCEPTION(
true,
Exceptions::RuntimeError,
"MueLu::ParameterListInterpreter():: invalid level id");
2197 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
2200 for (ParameterList::ConstIterator param = paramList.begin(); param != paramList.end(); ++param) {
2201 const std::string & paramName = paramList.name(param);
2202 const Teuchos::ParameterEntry & paramValue = paramList.entry(param);
2206 if (paramValue.isList()) {
2207 ParameterList paramList1 = Teuchos::getValue<ParameterList>(paramValue);
2208 if (paramList1.isParameter(
"factory")) {
2211 "MueLu::ParameterListInterpreter(): It seems that in the parameter lists for defining " << paramName <<
2212 " there is both a 'factory' and 'dependency for' parameter. This is not allowed. Please remove the 'dependency for' parameter.");
2214 factoryMapOut[paramName] = factFact_->BuildFactory(paramValue, factoryMapIn, factoryManagers);
2216 }
else if (paramList1.isParameter(
"dependency for")) {
2218 "MueLu::ParameterListInterpreter(): It seems that in the parameter lists for defining " << paramName <<
2219 " there is both a 'factory' and 'dependency for' parameter. This is not allowed.");
2221 std::string factoryName = paramList1.get<std::string>(
"dependency for");
2223 RCP<const FactoryBase> factbase = factoryMapIn.find(factoryName )->second;
2225 "MueLu::ParameterListInterpreter(): could not find factory " + factoryName +
" in factory map. Did you define it before?");
2227 RCP<const Factory> factoryconst = Teuchos::rcp_dynamic_cast<
const Factory>(factbase);
2228 RCP< Factory> factory = Teuchos::rcp_const_cast<
Factory>(factoryconst);
2232 for (ParameterList::ConstIterator vparam = validParamList->begin(); vparam != validParamList->end(); ++vparam) {
2233 const std::string& pName = validParamList->name(vparam);
2235 if (!paramList1.isParameter(pName)) {
2240 if (validParamList->isType< RCP<const FactoryBase> >(pName)) {
2242 RCP<const FactoryBase> generatingFact = factFact_->BuildFactory(paramList1.getEntry(pName), factoryMapIn, factoryManagers);
2243 factory->SetFactory(pName, generatingFact.create_weak());
2245 }
else if (validParamList->isType<RCP<const ParameterList> >(pName)) {
2246 if (pName ==
"ParameterList") {
2251 RCP<const ParameterList> subList = Teuchos::sublist(rcp(
new ParameterList(paramList1)), pName);
2252 factory->SetParameter(pName, ParameterEntry(subList));
2255 factory->SetParameter(pName, paramList1.getEntry(pName));
2259 }
else if (paramList1.isParameter(
"group")) {
2261 std::string groupType = paramList1.get<std::string>(
"group");
2263 "group must be of type \"FactoryManager\".");
2265 ParameterList groupList = paramList1;
2266 groupList.remove(
"group");
2268 bool setKokkosRefactor =
false;
2269 bool kokkosRefactor = useKokkos_;
2270 if (groupList.isParameter(
"use kokkos refactor")) {
2271 kokkosRefactor = groupList.get<
bool>(
"use kokkos refactor");
2272 groupList.remove(
"use kokkos refactor");
2273 setKokkosRefactor =
true;
2277 BuildFactoryMap(groupList, factoryMapIn, groupFactoryMap, factoryManagers);
2281 RCP<FactoryManager> m = rcp(
new FactoryManager(groupFactoryMap));
2282 if (setKokkosRefactor)
2283 m->SetKokkosRefactor(kokkosRefactor);
2284 factoryManagers[paramName] = m;
2287 this->GetOStream(
Warnings0) <<
"Could not interpret parameter list " << paramList1 << std::endl;
2289 "XML Parameter list must either be of type \"factory\" or of type \"group\".");
2293 factoryMapOut[paramName] = factFact_->BuildFactory(paramValue, factoryMapIn, factoryManagers);
2301 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
2304 Matrix& A =
dynamic_cast<Matrix&
>(Op);
2305 if (A.GetFixedBlockSize() != blockSize_)
2306 this->GetOStream(
Warnings0) <<
"Setting matrix block size to " << blockSize_ <<
" (value of the parameter in the list) " 2307 <<
"instead of " << A.GetFixedBlockSize() <<
" (provided matrix)." << std::endl
2308 <<
"You may want to check \"number of equations\" (or \"PDE equations\" for factory style list) parameter." << std::endl;
2310 A.SetFixedBlockSize(blockSize_, dofOffset_);
2312 #ifdef HAVE_MUELU_DEBUG 2313 MatrixUtils::checkLocalRowMapMatchesColMap(A);
2314 #endif // HAVE_MUELU_DEBUG 2316 }
catch (std::bad_cast& e) {
2317 this->GetOStream(
Warnings0) <<
"Skipping setting block size as the operator is not a matrix" << std::endl;
2321 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
2329 static bool compare(
const ParameterList& list1,
const ParameterList& list2) {
2332 for (ParameterList::ConstIterator it = list1.begin(); it != list1.end(); it++) {
2333 const std::string& name = it->first;
2334 const Teuchos::ParameterEntry& entry1 = it->second;
2336 const Teuchos::ParameterEntry *entry2 = list2.getEntryPtr(name);
2339 if (entry1.isList() && entry2->isList()) {
2340 compare(Teuchos::getValue<ParameterList>(entry1), Teuchos::getValue<ParameterList>(*entry2));
2343 if (entry1.getAny(
false) != entry2->getAny(
false))
2350 static inline bool areSame(
const ParameterList& list1,
const ParameterList& list2) {
2356 #define MUELU_PARAMETERLISTINTERPRETER_SHORT Important warning messages (one line)
Generic Smoother Factory for generating the smoothers of the MG hierarchy.
virtual RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
This class specifies the default factory that should generate some data on a Level if the data does n...
static void SetMueLuOFileStream(const std::string &filename)
Factory for determing the number of partitions for rebalancing.
#define MUELU_KOKKOS_FACTORY_NO_DECL(varName, oldFactory, newFactory)
Factory for generating coarse level map. Used by TentativePFactory.
std::string toString(const T &what)
Little helper function to convert non-string types to strings.
Factory for building transfer operators based on coarsening in polynomial degree, following the Intre...
Factory for building coarse grid matrices, when the matrix is of the form K+a*M. Useful when you want...
Factory that can generate other factories from.
void UpdateFactoryManager_CoarseSolvers(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps) const
static void SetDefaultVerbLevel(const VerbLevel defaultVerbLevel)
Set the default (global) verbosity level.
void UpdateFactoryManager_Smoothers(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps) const
Class that encapsulates external library smoothers.
size_t LevenshteinDistance(const char *s, size_t len_s, const char *t, size_t len_t)
static void DisableMultipleCheckGlobally()
#define MUELU_TEST_AND_SET_VAR(paramList, paramName, paramType, varName)
Factory for building permutation matrix that can be be used to shuffle data (matrices, vectors) among processes.
One-liner description of what is happening.
ParameterListInterpreter()
Empty constructor.
void BuildFactoryMap(const Teuchos::ParameterList ¶mList, const FactoryMap &factoryMapIn, FactoryMap &factoryMapOut, FactoryManagerMap &factoryManagers) const
Interpret "Factories" sublist.
void UpdateFactoryManager_SA(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps) const
virtual void SetupOperator(Operator &A) const
Setup Operator object.
void SetCycle(CycleType Cycle)
Supports VCYCLE and WCYCLE types.
Namespace for MueLu classes and methods.
#define MUELU_KOKKOS_FACTORY(varName, oldFactory, newFactory)
Interface to Zoltan library.This interface provides access to partitioning methods in Zoltan...
bool IsParamMuemexVariable(const std::string &name)
Factory for creating a graph base on a given matrix.
void UpdateFactoryManager_Repartition(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps, RCP< Factory > &nullSpaceFactory) const
static VerbLevel GetDefaultVerbLevel()
Get the default (global) verbosity level.
void UpdateFactoryManager_Aggregation_TentativeP(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps) const
std::map< std::string, RCP< FactoryManagerBase > > FactoryManagerMap
std::map< std::string, RCP< const FactoryBase > > FactoryMap
void UpdateFactoryManager_PG(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps) const
void UpdateFactoryManager_RAP(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps) const
Factory for building tentative prolongator.
#define TEST_MUTUALLY_EXCLUSIVE(arg1, arg2)
MsgType toVerbLevel(const std::string &verbLevelStr)
Class that encapsulates direct solvers. Autoselection of AmesosSmoother or Amesos2Smoother according ...
Factory for coarsening a graph with uncoupled aggregation.
virtual void SetupHierarchy(Hierarchy &H) const
Setup Hierarchy object.
Prolongator factory performing semi-coarsening.
Factory for building restriction operators using a prolongator factory.
void UpdateFactoryManager_PCoarsen(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps) const
Teuchos::RCP< MueLu::FacadeClassFactory< Scalar, LocalOrdinal, GlobalOrdinal, Node > > facadeFact_
FacadeClass factory.
#define MUELU_TEST_AND_SET_PARAM_2LIST(paramList, defaultList, paramName, paramType, listWrite)
static bool compare(const ParameterList &list1, const ParameterList &list2)
static CycleType GetDefaultCycle()
void UpdateFactoryManager(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps) const
MueLu::DefaultGlobalOrdinal GlobalOrdinal
void SetCycleStartLevel(int cycleStart)
Factory for interacting with Matlab.
Factory for interacting with Matlab.
Factory for building line detection information.
void UpdateFactoryManager_Nullspace(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps, RCP< Factory > &nullSpaceFactory) const
void Validate(const Teuchos::ParameterList ¶mList) const
void SetFactoryParameterList(const Teuchos::ParameterList ¶mList)
Factory interpreter stuff.
Factory to export aggregation info or visualize aggregates using VTK.
Prolongator factory which allows switching between two different prolongator strategies.
Interface to Zoltan2 library.This interface provides access to partitioning methods in Zoltan2...
void UpdateFactoryManager_SemiCoarsen(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps) const
Factory for building the constraint operator.
void UpdateFactoryManager_Emin(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps) const
Applies permutation to grid transfer operators.
void SetFactory(const std::string &varName, const RCP< const FactoryBase > &factory)
Set Factory.
const RCP< FactoryBase > GetFactoryNonConst(const std::string &varName)
Get factory associated with a particular data name (NONCONST version)
const RCP< const FactoryBase > GetFactory(const std::string &varName) const
Get factory associated with a particular data name.
Factory for generating a very special nullspace.
static void EnableTimerSync()
void SetupHierarchy(Hierarchy &H) const
Call the SetupHierarchy routine from the HiearchyManager object.
Factory for creating a graph base on a given matrix.
Class that encapsulates Matlab smoothers.
Partitioning within a node onlyThis interface provides partitioning within a node.
void UpdateFactoryManager_Coordinates(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps) const
void SetParameterList(const Teuchos::ParameterList ¶mList)
Set parameter list for Parameter list interpreter.
Class for transferring coordinates from a finer level to a coarser one.
Factory for building nonzero patterns for energy minimization.
static Teuchos::RCP< Teuchos::ParameterList > GetProblemSpecificList(std::string const &problemType)
Return default parameter settings for the specified problem type.
void UpdateFactoryManager_Matlab(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps) const
Class for transferring coordinates from a finer level to a coarser one.
Factory for building restriction operators.
Factory for building Energy Minimization prolongators.
void UpdateFactoryManager_Restriction(Teuchos::ParameterList ¶mList, const Teuchos::ParameterList &defaultList, FactoryManager &manager, int levelID, std::vector< keep_pair > &keeps) const
static bool areSame(const ParameterList &list1, const ParameterList &list2)
Helper functions to compare two paramter lists.
#define TEST_MUTUALLY_EXCLUSIVE_S(arg1, arg2)
static int GetDefaultCycleStartLevel()
void SetProlongatorScalingFactor(double scalingFactor)
Specify damping factor alpha such that x = x + alpha*P*c, where c is the coarse grid correction...
#define MUELU_TEST_PARAM_2LIST(paramList, defaultList, paramName, paramType, cmpValue)
Factory for building coarse matrices.
Exception throws to report errors in the internal logical of the program.
Factory for building filtered matrices using filtered graphs.
std::pair< std::string, const FactoryBase * > keep_pair
Description of what is happening (more verbose)
Factory for building coarse matrices.
Factory for building Petrov-Galerkin Smoothed Aggregation prolongators.
An factory which assigns each aggregate a quality estimate. Originally developed by Napov and Notay i...
Factory for building Smoothed Aggregation prolongators.
Factory for building uncoupled aggregates.
static Teuchos::RCP< const Teuchos::ParameterList > List()
Return a "master" list of all valid parameters and their default values.
#define MUELU_SET_VAR_2LIST(paramList, defaultList, paramName, paramType, varName)
void SetEasyParameterList(const Teuchos::ParameterList ¶mList)
Provides methods to build a multigrid hierarchy and apply multigrid cycles.
Factory for generating nullspace.
long ExtractNonSerializableData(const Teuchos::ParameterList &inList, Teuchos::ParameterList &serialList, Teuchos::ParameterList &nonSerialList)
Extract non-serializable data from level-specific sublists and move it to a separate parameter list...
Exception throws to report invalid user entry.
static const RCP< const NoFactory > getRCP()
Static Get() functions.