diff --git a/packages/muelu/doc/UsersGuide/masterList.xml b/packages/muelu/doc/UsersGuide/masterList.xml
index 48fbba6e50a1..cd1ba275f8f5 100644
--- a/packages/muelu/doc/UsersGuide/masterList.xml
+++ b/packages/muelu/doc/UsersGuide/masterList.xml
@@ -781,6 +781,24 @@
eigen-analysis: iterations
+
+ gmg: interpolation order
+ int
+ 1
+ Interpolation order used to interpolate values from coarse points to fine points. Possible values are 0 for piece-wise constant interpolation and 1 for piece-wise linear interpolation. This parameter is set to 1 by default.
+
+ parameter not existing in ML
+
+
+
+ gmg: build coarse coordinates
+ bool
+ true
+ If false, skip the calculation of coarse coordinates.
+ false
+ not supported by ML
+
+
transfer: params
\parameterlist
@@ -797,7 +815,6 @@
false
-
pcoarsen: schedule
string
diff --git a/packages/muelu/doc/UsersGuide/paramlist_hidden.tex b/packages/muelu/doc/UsersGuide/paramlist_hidden.tex
index e18f02201a1c..d9d3a4f770b9 100644
--- a/packages/muelu/doc/UsersGuide/paramlist_hidden.tex
+++ b/packages/muelu/doc/UsersGuide/paramlist_hidden.tex
@@ -154,6 +154,10 @@
\cbb{sa: eigenvalue estimate num iterations}{int}{10}{Number of power iterations to use to estimate max eigenvalue. This has effect only if estimation occurs during prolongator smoothing.}
+\cbb{gmg: interpolation order}{int}{1}{Interpolation order used to interpolate values from coarse points to fine points. Possible values are 0 for piece-wise constant interpolation and 1 for piece-wise linear interpolation. This parameter is set to 1 by default.}
+
+\cbb{gmg: build coarse coordinates}{bool}{true}{If false, skip the calculation of coarse coordinates.}
+
\cba{transfer: params}{\parameterlist}{Sublist of options for use by transfer.}
\cbb{pcoarsen: element}{string}{""}{Element class for IntrepidPCoarsening}
diff --git a/packages/muelu/src/CMakeLists.txt b/packages/muelu/src/CMakeLists.txt
index c627c046cd56..f35738136d96 100644
--- a/packages/muelu/src/CMakeLists.txt
+++ b/packages/muelu/src/CMakeLists.txt
@@ -472,3 +472,4 @@ TRIBITS_ADD_LIBRARY(
# touch CMakeLists.txt because a new file was created in Utils/ExplicitInstantiation of Utils/ForwardDeclaration
# touch CMakeLists.txt because a new file was created in Utils/ExplicitInstantiation of Utils/ForwardDeclaration
# touch CMakeLists.txt because a new file was created in Utils/ExplicitInstantiation of Utils/ForwardDeclaration
+# touch CMakeLists.txt because a new file was created in Utils/ExplicitInstantiation of Utils/ForwardDeclaration
diff --git a/packages/muelu/src/Graph/MueLu_Aggregates_decl.hpp b/packages/muelu/src/Graph/MueLu_Aggregates_decl.hpp
index f125e82d3ba2..6fd5f4ceeec3 100644
--- a/packages/muelu/src/Graph/MueLu_Aggregates_decl.hpp
+++ b/packages/muelu/src/Graph/MueLu_Aggregates_decl.hpp
@@ -136,7 +136,7 @@ namespace MueLu {
This has to be done by the aggregation factory.
*/
- RCP & GetIndexManager() { return geoData_; }
+ RCP& GetIndexManager() { return geoData_; }
/*! @brief Get the index manager used by structured aggregation algorithms.
diff --git a/packages/muelu/src/Graph/StructuredAggregation/MueLu_AggregationStructuredAlgorithm_decl.hpp b/packages/muelu/src/Graph/StructuredAggregation/MueLu_AggregationStructuredAlgorithm_decl.hpp
index 46f6b1e17c4c..b17d01ae2173 100644
--- a/packages/muelu/src/Graph/StructuredAggregation/MueLu_AggregationStructuredAlgorithm_decl.hpp
+++ b/packages/muelu/src/Graph/StructuredAggregation/MueLu_AggregationStructuredAlgorithm_decl.hpp
@@ -100,17 +100,28 @@ namespace MueLu {
/*! @brief Local aggregation. */
- void BuildAggregates(const ParameterList& params, const GraphBase& graph, Aggregates& aggregates,
- std::vector& aggStat, LO& numNonAggregatedNodes) const;
+ void BuildAggregates(const Teuchos::ParameterList& params, const GraphBase& graph,
+ Aggregates& aggregates, std::vector& aggStat,
+ LO& numNonAggregatedNodes) const;
/*! @brief Local aggregation. */
- void BuildAggregates(const ParameterList& params, const GraphBase& graph,
- RCP& geoData, RCP& myGraph) const;
+ void BuildAggregates(const GraphBase& graph, RCP& geoData,
+ RCP& myGraph) const;
//@}
std::string description() const { return "Aggretation: structured algorithm"; }
+ private:
+
+ void ComputeGraphDataConstant(const GraphBase& graph, RCP& geoData,
+ const int numInterpolationPoints, ArrayRCP& nnzOnRow,
+ Array& rowPtr, Array& colIndex) const;
+
+ void ComputeGraphDataLinear(const GraphBase& graph, RCP& geoData,
+ const int numInterpolationPoints, ArrayRCP& nnzOnRow,
+ Array& rowPtr, Array& colIndex) const;
+
};
} //namespace MueLu
diff --git a/packages/muelu/src/Graph/StructuredAggregation/MueLu_AggregationStructuredAlgorithm_def.hpp b/packages/muelu/src/Graph/StructuredAggregation/MueLu_AggregationStructuredAlgorithm_def.hpp
index b09d32806f6d..5fe919d11f6c 100644
--- a/packages/muelu/src/Graph/StructuredAggregation/MueLu_AggregationStructuredAlgorithm_def.hpp
+++ b/packages/muelu/src/Graph/StructuredAggregation/MueLu_AggregationStructuredAlgorithm_def.hpp
@@ -67,7 +67,8 @@ namespace MueLu {
template
void AggregationStructuredAlgorithm::
- BuildAggregates(const ParameterList& params, const GraphBase& graph, Aggregates& aggregates, std::vector& aggStat,
+ BuildAggregates(const Teuchos::ParameterList& params, const GraphBase& graph,
+ Aggregates& aggregates, std::vector& aggStat,
LO& numNonAggregatedNodes) const {
Monitor m(*this, "BuildAggregates");
@@ -79,22 +80,21 @@ namespace MueLu {
out = Teuchos::getFancyOStream(rcp(new Teuchos::oblackholestream()));
}
- std::string coupling = params.get("aggregation: coupling");
- const bool coupled = (coupling == "coupled" ? true : false);
+ RCP geoData = aggregates.GetIndexManager();
+ const bool coupled = geoData->isAggregationCoupled();
ArrayRCP vertex2AggId = aggregates.GetVertex2AggId()->getDataNonConst(0);
ArrayRCP procWinner = aggregates.GetProcWinner() ->getDataNonConst(0);
- RCP geoData = aggregates.GetIndexManager();
- Array ghostedCoarseNodeCoarseLIDs;
- Array ghostedCoarseNodeCoarsePIDs;
- Array ghostedCoarseNodeCoarseGIDs;
+ Array ghostedCoarseNodeCoarseLIDs;
+ Array ghostedCoarseNodeCoarsePIDs;
+ Array ghostedCoarseNodeCoarseGIDs;
*out << "Extract data for ghosted nodes" << std::endl;
geoData->getGhostedNodesData(graph.GetDomainMap(), ghostedCoarseNodeCoarseLIDs,
ghostedCoarseNodeCoarsePIDs, ghostedCoarseNodeCoarseGIDs);
- LO rem;
+ LO rem, rate;
Array ghostedIdx(3), coarseIdx(3);
- LO ghostedCoarseNodeCoarseLID, aggId, rate;
+ LO ghostedCoarseNodeCoarseLID, aggId;
*out << "Loop over fine nodes and assign them to an aggregate and a rank" << std::endl;
for(LO nodeIdx = 0; nodeIdx < geoData->getNumLocalFineNodes(); ++nodeIdx) {
// Compute coarse ID associated with fine LID
@@ -109,9 +109,9 @@ namespace MueLu {
} else {
rate = geoData->getCoarseningEndRate(dim);
}
- if(rem > (rate / 2)) { ++coarseIdx[dim]; }
+ if(rem > (rate / 2)) {++coarseIdx[dim];}
if(coupled && (geoData->getStartGhostedCoarseNode(dim)*geoData->getCoarseningRate(dim)
- > geoData->getStartIndex(dim))) { --coarseIdx[dim]; }
+ > geoData->getStartIndex(dim))) {--coarseIdx[dim];}
}
geoData->getCoarseNodeGhostedLID(coarseIdx[0], coarseIdx[1], coarseIdx[2],
@@ -124,14 +124,12 @@ namespace MueLu {
--numNonAggregatedNodes;
} // Loop over fine points
-
- }
+ } // BuildAggregates()
template
void AggregationStructuredAlgorithm::
- BuildAggregates(const ParameterList& params, const GraphBase& graph, RCP& geoData,
- RCP& myGraph) const {
+ BuildAggregates(const GraphBase& graph, RCP& geoData, RCP& myGraph) const {
Monitor m(*this, "BuildAggregates");
RCP out;
@@ -142,8 +140,7 @@ namespace MueLu {
out = Teuchos::getFancyOStream(rcp(new Teuchos::oblackholestream()));
}
- std::string coupling = params.get("aggregation: coupling");
- const bool coupled = (coupling == "coupled" ? true : false);
+ const bool coupled = geoData->isAggregationCoupled();
Array ghostedCoarseNodeCoarseLIDs;
Array ghostedCoarseNodeCoarsePIDs;
Array ghostedCoarseNodeCoarseGIDs;
@@ -152,24 +149,138 @@ namespace MueLu {
geoData->getGhostedNodesData(graph.GetDomainMap(), ghostedCoarseNodeCoarseLIDs,
ghostedCoarseNodeCoarsePIDs, ghostedCoarseNodeCoarseGIDs);
- int numDimensions = geoData->getNumDimensions();
- Array ghostedIdx(3,0);
- Array coarseIdx(3,0);
- Array ijkRem(3,0);
-
- LO ghostedCoarseNodeCoarseLID, aggId, rate;
- GO myGID;
-
- bool allCoarse;
- int numInterpolationPoints = 1 << numDimensions; // Using bit logic to compute 2^numDimensions
+ // Compute the number of coarse points needed to interpolate quantities to a fine point
+ int numInterpolationPoints = 0;
+ if(geoData->getInterpolationOrder() == 0) {
+ numInterpolationPoints = 1;
+ } else if(geoData->getInterpolationOrder() == 1) {
+ // Compute 2^numDimensions using bit logic to avoid round-off errors
+ numInterpolationPoints = 1 << geoData->getNumDimensions();
+ }
*out << "numInterpolationPoints=" << numInterpolationPoints << std::endl;
- Array isCoarse(numDimensions);
- Array colIndex(geoData->getNumLocalFineNodes()*numInterpolationPoints);// also resize...
- Array rowPtr(geoData->getNumLocalFineNodes()+1);
+
+ Array colIndex( geoData->getNumLocalCoarseNodes() + numInterpolationPoints*
+ (geoData->getNumLocalFineNodes() - geoData->getNumLocalCoarseNodes()) );
+ Array rowPtr(geoData->getNumLocalFineNodes()+1);
rowPtr[0] = 0;
ArrayRCP nnzOnRow(geoData->getNumLocalFineNodes());
- *out << "Loop over fine nodes and assign them to an aggregate and a rank" << std::endl;
+ *out << "Compute prolongatorGraph data" << std::endl;
+ if(geoData->getInterpolationOrder() == 0) {
+ ComputeGraphDataConstant(graph, geoData, numInterpolationPoints, nnzOnRow, rowPtr, colIndex);
+ } else if(geoData->getInterpolationOrder() == 1) {
+ ComputeGraphDataLinear(graph, geoData, numInterpolationPoints, nnzOnRow, rowPtr, colIndex);
+ }
+
+ //Redo colMap using colData to specify the coarse nodes and ghost nodes needed
+ RCP