Skip to content

Commit

Permalink
Merge 'trilinos/Trilinos:develop' (fe80646) into 'EM-Plasma/Trilinos:…
Browse files Browse the repository at this point in the history
…develop' (076963d).

* trilinos/develop:
  Zoltan2 tpetra dep (trilinos#9287)
  MueLu SmootherFactory: Fix wrong ""PostSmoother data""
  MueLu Maxwell test: Allow reading in map files
  MueLu RefMaxwell: Label MultiVectors
  MueLu Maxwell test: Fix reuse
  MueLu RefMaxwell: Overlap imports into subproblems
  Xpetra: Expose Tpetra's {begin,end}{Import,Export} for MVs
  RefMaxwell: Allow reuse of fine level smoother
  • Loading branch information
EMPIRE Jenkins committed Jun 30, 2021
2 parents 076963d + fe80646 commit 3fe9124
Show file tree
Hide file tree
Showing 29 changed files with 367 additions and 143 deletions.
1 change: 1 addition & 0 deletions packages/muelu/adapters/xpetra/MueLu_RefMaxwell_decl.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -421,6 +421,7 @@ namespace MueLu {
//! Two hierarchies: one for the coarse (1,1)-block, another for the (2,2)-block
Teuchos::RCP<Hierarchy> HierarchyH_, Hierarchy22_;
Teuchos::RCP<SmootherBase> PreSmoother_, PostSmoother_;
Teuchos::RCP<SmootherPrototype> PreSmootherData_, PostSmootherData_;
#if defined(MUELU_REFMAXWELL_CAN_USE_HIPTMAIR)
Teuchos::RCP<Ifpack2::Preconditioner<Scalar,LocalOrdinal,GlobalOrdinal,Node> > hiptmairPreSmoother_, hiptmairPostSmoother_;
#endif
Expand Down
103 changes: 71 additions & 32 deletions packages/muelu/adapters/xpetra/MueLu_RefMaxwell_def.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -1324,6 +1324,15 @@ namespace MueLu {
throw(Xpetra::Exceptions::RuntimeError("MueLu must be compiled with Ifpack2 for Hiptmair smoothing."));
#endif // defined(MUELU_REFMAXWELL_CAN_USE_HIPTMAIR)
} else {

Level level;
RCP<MueLu::FactoryManagerBase> factoryHandler = rcp(new FactoryManager());
level.SetFactoryManager(factoryHandler);
level.SetLevelID(0);
level.setObjectLabel("RefMaxwell (1,1)");
level.Set("A",SM_Matrix_);
level.setlib(SM_Matrix_->getDomainMap()->lib());

if (parameterList_.isType<std::string>("smoother: pre type") && parameterList_.isType<std::string>("smoother: post type")) {
std::string preSmootherType = parameterList_.get<std::string>("smoother: pre type");
std::string postSmootherType = parameterList_.get<std::string>("smoother: post type");
Expand All @@ -1334,42 +1343,49 @@ namespace MueLu {
if (parameterList_.isSublist("smoother: post params"))
postSmootherList = parameterList_.sublist("smoother: post params");

Level level;
RCP<MueLu::FactoryManagerBase> factoryHandler = rcp(new FactoryManager());
level.SetFactoryManager(factoryHandler);
level.SetLevelID(0);
level.setObjectLabel("RefMaxwell (1,1)");
level.Set("A",SM_Matrix_);
level.setlib(SM_Matrix_->getDomainMap()->lib());

RCP<SmootherPrototype> preSmootherPrototype = rcp(new TrilinosSmoother(preSmootherType, preSmootherList));
RCP<SmootherFactory> preSmootherFact = rcp(new SmootherFactory(preSmootherPrototype));

RCP<SmootherPrototype> postSmootherPrototype = rcp(new TrilinosSmoother(postSmootherType, postSmootherList));
RCP<SmootherFactory> postSmootherFact = rcp(new SmootherFactory(postSmootherPrototype));

level.Request("PreSmoother",preSmootherFact.get());
preSmootherFact->Build(level);
PreSmoother_ = level.Get<RCP<SmootherBase> >("PreSmoother",preSmootherFact.get());
RCP<SmootherFactory> smootherFact = rcp(new SmootherFactory(preSmootherPrototype, postSmootherPrototype));

level.Request("PostSmoother",postSmootherFact.get());
postSmootherFact->Build(level);
PostSmoother_ = level.Get<RCP<SmootherBase> >("PostSmoother",postSmootherFact.get());
level.Request("PreSmoother",smootherFact.get());
level.Request("PostSmoother",smootherFact.get());
if (enable_reuse_) {
ParameterList smootherFactoryParams;
smootherFactoryParams.set("keep smoother data", true);
smootherFact->SetParameterList(smootherFactoryParams);
level.Request("PreSmoother data", smootherFact.get());
level.Request("PostSmoother data", smootherFact.get());
if (!PreSmootherData_.is_null())
level.Set("PreSmoother data", PreSmootherData_, smootherFact.get());
if (!PostSmootherData_.is_null())
level.Set("PostSmoother data", PostSmootherData_, smootherFact.get());
}
smootherFact->Build(level);
PreSmoother_ = level.Get<RCP<SmootherBase> >("PreSmoother",smootherFact.get());
PostSmoother_ = level.Get<RCP<SmootherBase> >("PostSmoother",smootherFact.get());
if (enable_reuse_) {
PreSmootherData_ = level.Get<RCP<SmootherPrototype> >("PreSmoother data",smootherFact.get());
PostSmootherData_ = level.Get<RCP<SmootherPrototype> >("PostSmoother data",smootherFact.get());
}
} else {
std::string smootherType = parameterList_.get<std::string>("smoother: type", "CHEBYSHEV");
Level level;
RCP<MueLu::FactoryManagerBase> factoryHandler = rcp(new FactoryManager());
level.SetFactoryManager(factoryHandler);
level.SetLevelID(0);
level.setObjectLabel("RefMaxwell (1,1)");
level.Set("A",SM_Matrix_);
level.setlib(SM_Matrix_->getDomainMap()->lib());

RCP<SmootherPrototype> smootherPrototype = rcp(new TrilinosSmoother(smootherType, smootherList_));
RCP<SmootherFactory> SmootherFact = rcp(new SmootherFactory(smootherPrototype));
level.Request("PreSmoother",SmootherFact.get());
SmootherFact->Build(level);
PreSmoother_ = level.Get<RCP<SmootherBase> >("PreSmoother",SmootherFact.get());
RCP<SmootherFactory> smootherFact = rcp(new SmootherFactory(smootherPrototype));
level.Request("PreSmoother",smootherFact.get());
if (enable_reuse_) {
ParameterList smootherFactoryParams;
smootherFactoryParams.set("keep smoother data", true);
smootherFact->SetParameterList(smootherFactoryParams);
level.Request("PreSmoother data", smootherFact.get());
if (!PreSmootherData_.is_null())
level.Set("PreSmoother data", PreSmootherData_, smootherFact.get());
}
smootherFact->Build(level);
PreSmoother_ = level.Get<RCP<SmootherBase> >("PreSmoother",smootherFact.get());
PostSmoother_ = PreSmoother_;
if (enable_reuse_)
PreSmootherData_ = level.Get<RCP<SmootherPrototype> >("PreSmoother data",smootherFact.get());
}
useHiptmairSmoothing_ = false;
}
Expand All @@ -1385,43 +1401,56 @@ namespace MueLu {
P11res_ = MultiVectorFactory::Build(R11_->getRangeMap(), numVectors);
else
P11res_ = MultiVectorFactory::Build(P11_->getDomainMap(), numVectors);
if (D0_T_R11_colMapsMatch_)
P11res_->setObjectLabel("P11res");
if (D0_T_R11_colMapsMatch_) {
D0TR11Tmp_ = MultiVectorFactory::Build(R11_->getColMap(), numVectors);
D0TR11Tmp_->setObjectLabel("D0TR11Tmp");
}
if (!ImporterH_.is_null()) {
P11resTmp_ = MultiVectorFactory::Build(ImporterH_->getTargetMap(), numVectors);
P11resTmp_->setObjectLabel("P11resTmp");
P11x_ = MultiVectorFactory::Build(ImporterH_->getTargetMap(), numVectors);
} else
P11x_ = MultiVectorFactory::Build(P11_->getDomainMap(), numVectors);
P11x_->setObjectLabel("P11x");
if (!D0_T_Matrix_.is_null())
D0res_ = MultiVectorFactory::Build(D0_T_Matrix_->getRangeMap(), numVectors);
else
D0res_ = MultiVectorFactory::Build(D0_Matrix_->getDomainMap(), numVectors);
D0res_->setObjectLabel("D0res");
if (!Importer22_.is_null()) {
D0resTmp_ = MultiVectorFactory::Build(Importer22_->getTargetMap(), numVectors);
D0resTmp_->setObjectLabel("D0resTmp");
D0x_ = MultiVectorFactory::Build(Importer22_->getTargetMap(), numVectors);
} else
D0x_ = MultiVectorFactory::Build(D0_Matrix_->getDomainMap(), numVectors);
D0x_->setObjectLabel("D0x");
if (!AH_.is_null()) {
if (!ImporterH_.is_null() && !implicitTranspose_)
P11resSubComm_ = MultiVectorFactory::Build(P11resTmp_, Teuchos::View);
else
P11resSubComm_ = MultiVectorFactory::Build(P11res_, Teuchos::View);
P11resSubComm_->replaceMap(AH_->getRangeMap());
P11resSubComm_->setObjectLabel("P11resSubComm");

P11xSubComm_ = MultiVectorFactory::Build(P11x_, Teuchos::View);
P11xSubComm_->replaceMap(AH_->getDomainMap());
P11xSubComm_->setObjectLabel("P11xSubComm");
}
if (!A22_.is_null()) {
if (!Importer22_.is_null() && !implicitTranspose_)
D0resSubComm_ = MultiVectorFactory::Build(D0resTmp_, Teuchos::View);
else
D0resSubComm_ = MultiVectorFactory::Build(D0res_, Teuchos::View);
D0resSubComm_->replaceMap(A22_->getRangeMap());
D0resSubComm_->setObjectLabel("D0resSubComm");

D0xSubComm_ = MultiVectorFactory::Build(D0x_, Teuchos::View);
D0xSubComm_->replaceMap(A22_->getDomainMap());
D0xSubComm_->setObjectLabel("D0xSubComm");
}
residual_ = MultiVectorFactory::Build(SM_Matrix_->getDomainMap(), numVectors);
residual_->setObjectLabel("residual");
}


Expand Down Expand Up @@ -2446,15 +2475,18 @@ namespace MueLu {

if (!ImporterH_.is_null() && !implicitTranspose_) {
RCP<Teuchos::TimeMonitor> tmH = getTimer("MueLu RefMaxwell: import coarse (1,1)");
P11resTmp_->doImport(*P11res_, *ImporterH_, Xpetra::INSERT);
P11resTmp_->beginImport(*P11res_, *ImporterH_, Xpetra::INSERT);
}
if (!allNodesBoundary_ && !Importer22_.is_null() && !implicitTranspose_) {
RCP<Teuchos::TimeMonitor> tm22 = getTimer("MueLu RefMaxwell: import (2,2)");
D0resTmp_->doImport(*D0res_, *Importer22_, Xpetra::INSERT);
D0resTmp_->beginImport(*D0res_, *Importer22_, Xpetra::INSERT);
}

// iterate on coarse (1, 1) block
if (!AH_.is_null()) {
if (!ImporterH_.is_null() && !implicitTranspose_)
P11resTmp_->endImport(*P11res_, *ImporterH_, Xpetra::INSERT);

RCP<Teuchos::TimeMonitor> tmH = getTimer("MueLu RefMaxwell: solve coarse (1,1)", AH_->getRowMap()->getComm());

#if defined(HAVE_MUELU_STRATIMIKOS) && defined(HAVE_MUELU_THYRA)
Expand All @@ -2473,6 +2505,9 @@ namespace MueLu {

// iterate on (2, 2) block
if (!A22_.is_null()) {
if (!allNodesBoundary_ && !Importer22_.is_null() && !implicitTranspose_)
D0resTmp_->endImport(*D0res_, *Importer22_, Xpetra::INSERT);

RCP<Teuchos::TimeMonitor> tm22 = getTimer("MueLu RefMaxwell: solve (2,2)", A22_->getRowMap()->getComm());

#if defined(HAVE_MUELU_STRATIMIKOS) && defined(HAVE_MUELU_THYRA)
Expand All @@ -2489,6 +2524,10 @@ namespace MueLu {
Hierarchy22_->Iterate(*D0resSubComm_, *D0xSubComm_, numIters22_, true);
}

if (AH_.is_null() && !ImporterH_.is_null() && !implicitTranspose_)
P11resTmp_->endImport(*P11res_, *ImporterH_, Xpetra::INSERT);
if (A22_.is_null() && !allNodesBoundary_ && !Importer22_.is_null() && !implicitTranspose_)
D0resTmp_->endImport(*D0res_, *Importer22_, Xpetra::INSERT);
}

if (fuseProlongationAndUpdate_) {
Expand Down
4 changes: 2 additions & 2 deletions packages/muelu/src/Smoothers/MueLu_SmootherFactory_def.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -126,7 +126,7 @@ namespace MueLu {
void SmootherFactory<Scalar, LocalOrdinal, GlobalOrdinal, Node>::BuildSmoother(Level& currentLevel, PreOrPost const preOrPost) const {
// SmootherFactory is quite tricky because of the fact that one of the smoother prototypes may be zero.
// The challenge is that we have no way of knowing how user uses this factory. For instance, lets say
// user wants to use s1 prototype as a presmoother, and s2 as a postsmoother. He could do:
// user wants to use s1 prototype as a presmoother, and s2 as a postsmoother. They could do:
// (a) create SmootherFactory(s1, s2), or
// (b) create SmootherFactory(s1, null) and SmootherFactory(null, s2)
// It may also happen that somewhere somebody set presmoother factory = postsmoother factory = (a)
Expand Down Expand Up @@ -227,7 +227,7 @@ namespace MueLu {
currentLevel.Set<RCP<SmootherBase> >("PostSmoother", postSmoother, this);

if (pL.get<bool>("keep smoother data"))
Set(currentLevel, "PostSmoother data", preSmoother);
Set(currentLevel, "PostSmoother data", postSmoother);
}

ParameterList& paramList = const_cast<ParameterList&>(this->GetParameterList());
Expand Down
41 changes: 36 additions & 5 deletions packages/muelu/test/maxwell/Maxwell3D.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -381,6 +381,13 @@ bool SetupSolveWrappers<Scalar,LocalOrdinal,GlobalOrdinal,Node>::SetupSolve(std:
if (reuse) {
TEUCHOS_ASSERT(precType == "MueLu-RefMaxwell");
for (int solveno = 0; solveno<2; solveno++) {
// SM_Matrix->resumeFill();
// SM_Matrix->fillComplete();
if (X0.is_null())
X->putScalar(Teuchos::ScalarTraits<Scalar>::zero());
else
X = X0;
problem -> setProblem( X, B );
Teuchos::rcp_dynamic_cast<MueLu::RefMaxwell<SC,LO,GO,NO> >(preconditioner)->resetMatrix(SM_Matrix);
Belos::ReturnType status = solver -> solve();
int iters = solver -> getNumIters();
Expand Down Expand Up @@ -516,6 +523,13 @@ bool SetupSolveWrappers<double,LocalOrdinal,GlobalOrdinal,Node>::SetupSolve(std:
if (reuse) {
TEUCHOS_ASSERT(precType == "MueLu-RefMaxwell");
for (int solveno = 0; solveno<2; solveno++) {
// SM_Matrix->resumeFill();
// SM_Matrix->fillComplete();
if (X0.is_null())
X->putScalar(Teuchos::ScalarTraits<Scalar>::zero());
else
X = X0;
problem -> setProblem( X, B );
Teuchos::rcp_dynamic_cast<MueLu::RefMaxwell<SC,LO,GO,NO> >(preconditioner)->resetMatrix(SM_Matrix);
Belos::ReturnType status = solver -> solve();
int iters = solver -> getNumIters();
Expand Down Expand Up @@ -730,12 +744,29 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib lib, int arg
// Read matrices in from files
RCP<Matrix> D0_Matrix, SM_Matrix, M1_Matrix, Ms_Matrix, M0inv_Matrix, Kn_Matrix;

// gradient matrix
D0_Matrix = Xpetra::IO<SC, LO, GO, NO>::Read(D0_file, lib, comm);

// maps for nodal and edge matrices
RCP<const Map> node_map = D0_Matrix->getDomainMap();
RCP<const Map> edge_map = D0_Matrix->getRangeMap();
RCP<const Map> node_map;
RCP<const Map> edge_map;

// gradient matrix
try {
std::string base = D0_file.substr(0, D0_file.find_last_of('/')+1);
std::string D0_filename = D0_file.substr(D0_file.find_last_of('/')+1, std::string::npos);
std::string edgeMap_file = base + "rowmap_" + D0_filename;
std::string nodeMap_file = base + "domainmap_" + D0_filename;
std::string colMap_file = base + "colmap_" + D0_filename;
node_map = Xpetra::IO<SC, LO, GO, NO>::ReadMap(nodeMap_file, lib, comm);
edge_map = Xpetra::IO<SC, LO, GO, NO>::ReadMap(edgeMap_file, lib, comm);
RCP<const Map> colMap;
if (comm->getSize() > 1)
colMap = Xpetra::IO<SC, LO, GO, NO>::ReadMap(colMap_file, lib, comm);
D0_Matrix = Xpetra::IO<SC, LO, GO, NO>::Read(D0_file, edge_map, colMap, node_map, edge_map);
} catch (const std::exception& e) {
// *out << "Skipping D0 maps, because: " << e.what() << std::endl;
D0_Matrix = Xpetra::IO<SC, LO, GO, NO>::Read(D0_file, lib, comm);
node_map = D0_Matrix->getDomainMap();
edge_map = D0_Matrix->getRangeMap();
}

// build stiffness plus mass matrix (SM_Matrix)
if (SM_file == "") {
Expand Down
26 changes: 25 additions & 1 deletion packages/xpetra/src/DistObject/Xpetra_DistObject.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -81,14 +81,38 @@ namespace Xpetra {
//! Import data into this object using an Import object ("forward mode").
virtual void doImport(const DistObject< Packet, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)= 0;

//! Import data into this object using an Import object ("forward mode").
virtual void beginImport(const DistObject< Packet, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM) { doImport(source, importer, CM); };

//! Import data into this object using an Import object ("forward mode").
virtual void endImport(const DistObject< Packet, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM) { };

//! Export data into this object using an Export object ("forward mode").
virtual void doExport(const DistObject< Packet, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)= 0;

//! Export data into this object using an Export object ("forward mode").
virtual void beginExport(const DistObject< Packet, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM) { doExport(source, exporter, CM); };

//! Export data into this object using an Export object ("forward mode").
virtual void endExport(const DistObject< Packet, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM) { };

//! Import data into this object using an Export object ("reverse mode").
virtual void doImport(const DistObject< Packet, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)= 0;

//! Import data into this object using an Export object ("reverse mode").
virtual void beginImport(const DistObject< Packet, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM) { doImport(source, exporter, CM); };

//! Import data into this object using an Export object ("reverse mode").
virtual void endImport(const DistObject< Packet, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM) { };

//! Export data into this object using an Import object ("reverse mode").
virtual void doExport(const DistObject< Packet, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)=0;

//! Export data into this object using an Import object ("reverse mode").
virtual void beginExport(const DistObject< Packet, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM) { doExport(source, importer, CM); };

//! Export data into this object using an Import object ("reverse mode").
virtual void doExport(const DistObject< Packet, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)= 0;
virtual void endExport(const DistObject< Packet, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM) { };

//@}

Expand Down
16 changes: 16 additions & 0 deletions packages/xpetra/src/MultiVector/Xpetra_TpetraMultiVector_decl.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -255,12 +255,28 @@ namespace Xpetra {

void doImport(const DistObject< Scalar, LocalOrdinal,GlobalOrdinal,Node> &source, const Import<LocalOrdinal,GlobalOrdinal,Node> &importer, CombineMode CM);

void beginImport(const DistObject< Scalar, LocalOrdinal,GlobalOrdinal,Node> &source, const Import<LocalOrdinal,GlobalOrdinal,Node> &importer, CombineMode CM);

void endImport(const DistObject< Scalar, LocalOrdinal,GlobalOrdinal,Node> &source, const Import<LocalOrdinal,GlobalOrdinal,Node> &importer, CombineMode CM);

void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import<LocalOrdinal,GlobalOrdinal,Node>& importer, CombineMode CM);

void beginExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import<LocalOrdinal,GlobalOrdinal,Node>& importer, CombineMode CM);

void endExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import<LocalOrdinal,GlobalOrdinal,Node>& importer, CombineMode CM);

void doImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Export<LocalOrdinal,GlobalOrdinal,Node>& exporter, CombineMode CM);

void beginImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Export<LocalOrdinal,GlobalOrdinal,Node>& exporter, CombineMode CM);

void endImport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &source, const Export<LocalOrdinal,GlobalOrdinal,Node>& exporter, CombineMode CM);

void doExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export<LocalOrdinal,GlobalOrdinal,Node>& exporter, CombineMode CM);

void beginExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export<LocalOrdinal,GlobalOrdinal,Node>& exporter, CombineMode CM);

void endExport(const DistObject< Scalar, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export<LocalOrdinal,GlobalOrdinal,Node>& exporter, CombineMode CM);

void replaceMap(const RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& map);

//@}
Expand Down
Loading

0 comments on commit 3fe9124

Please sign in to comment.