46 #ifndef MUELU_COORDINATESTRANSFER_FACTORY_DEF_HPP 47 #define MUELU_COORDINATESTRANSFER_FACTORY_DEF_HPP 49 #include "Xpetra_ImportFactory.hpp" 50 #include "Xpetra_MultiVectorFactory.hpp" 51 #include "Xpetra_MapFactory.hpp" 52 #include "Xpetra_IO.hpp" 54 #include "MueLu_CoarseMapFactory.hpp" 55 #include "MueLu_Aggregates.hpp" 64 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
66 RCP<ParameterList> validParamList = rcp(
new ParameterList());
68 validParamList->set<RCP<const FactoryBase> >(
"Coordinates", Teuchos::null,
"Factory for coordinates generation");
69 validParamList->set<RCP<const FactoryBase> >(
"Aggregates", Teuchos::null,
"Factory for coordinates generation");
70 validParamList->set<RCP<const FactoryBase> >(
"CoarseMap", Teuchos::null,
"Generating factory of the coarse map");
71 validParamList->set<
bool> (
"structured aggregation",
false,
"Flag specifying that the geometric data is transferred for StructuredAggregationFactory");
72 validParamList->set<
bool> (
"aggregation coupled",
false,
"Flag specifying if the aggregation algorithm was used in coupled mode.");
73 validParamList->set<
bool> (
"Geometric",
false,
"Flag specifying that the coordinates are transferred for GeneralGeometricPFactory");
74 validParamList->set<RCP<const FactoryBase> >(
"coarseCoordinates", Teuchos::null,
"Factory for coarse coordinates generation");
75 validParamList->set<RCP<const FactoryBase> >(
"gCoarseNodesPerDim", Teuchos::null,
"Factory providing the global number of nodes per spatial dimensions of the mesh");
76 validParamList->set<RCP<const FactoryBase> >(
"lCoarseNodesPerDim", Teuchos::null,
"Factory providing the local number of nodes per spatial dimensions of the mesh");
77 validParamList->set<RCP<const FactoryBase> >(
"numDimensions" , Teuchos::null,
"Factory providing the number of spatial dimensions of the mesh");
78 validParamList->set<
int> (
"write start", -1,
"first level at which coordinates should be written to file");
79 validParamList->set<
int> (
"write end", -1,
"last level at which coordinates should be written to file");
80 validParamList->set<
bool> (
"hybrid aggregation",
false,
"Flag specifying that hybrid aggregation data is transfered for HybridAggregationFactory");
81 validParamList->set<RCP<const FactoryBase> >(
"aggregationRegionTypeCoarse", Teuchos::null,
"Factory indicating what aggregation type is to be used on the coarse level of the region");
82 validParamList->set<
bool> (
"interface aggregation",
false,
"Flag specifying that interface aggregation data is transfered for HybridAggregationFactory");
83 validParamList->set<RCP<const FactoryBase> >(
"coarseInterfacesDimensions", Teuchos::null,
"Factory providing coarseInterfacesDimensions");
84 validParamList->set<RCP<const FactoryBase> >(
"nodeOnCoarseInterface", Teuchos::null,
"Factory providing nodeOnCoarseInterface");
87 return validParamList;
90 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
92 static bool isAvailableCoords =
false;
94 const ParameterList& pL = GetParameterList();
95 if(pL.get<
bool>(
"structured aggregation") ==
true) {
96 if(pL.get<
bool>(
"aggregation coupled") ==
true) {
97 Input(fineLevel,
"gCoarseNodesPerDim");
99 Input(fineLevel,
"lCoarseNodesPerDim");
100 Input(fineLevel,
"numDimensions");
101 }
else if(pL.get<
bool>(
"Geometric") ==
true) {
102 Input(coarseLevel,
"coarseCoordinates");
103 Input(coarseLevel,
"gCoarseNodesPerDim");
104 Input(coarseLevel,
"lCoarseNodesPerDim");
105 }
else if(pL.get<
bool>(
"hybrid aggregation") ==
true) {
106 Input(fineLevel,
"aggregationRegionTypeCoarse");
107 Input(fineLevel,
"lCoarseNodesPerDim");
108 Input(fineLevel,
"numDimensions");
109 if(pL.get<
bool>(
"interface aggregation") ==
true) {
110 Input(fineLevel,
"coarseInterfacesDimensions");
111 Input(fineLevel,
"nodeOnCoarseInterface");
115 isAvailableCoords = coarseLevel.
IsAvailable(
"Coordinates",
this);
117 if (isAvailableCoords ==
false) {
118 Input(fineLevel,
"Coordinates");
119 Input(fineLevel,
"Aggregates");
120 Input(fineLevel,
"CoarseMap");
125 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
129 using xdMV = Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::coordinateType,LO,GO,NO>;
131 GetOStream(
Runtime0) <<
"Transferring coordinates" << std::endl;
134 RCP<xdMV> coarseCoords;
135 RCP<xdMV> fineCoords;
136 Array<GO> gCoarseNodesPerDir;
137 Array<LO> lCoarseNodesPerDir;
139 const ParameterList& pL = GetParameterList();
141 if(pL.get<
bool>(
"hybrid aggregation") ==
true) {
142 std::string regionType = Get<std::string>(fineLevel,
"aggregationRegionTypeCoarse");
143 numDimensions = Get<int>(fineLevel,
"numDimensions");
144 lCoarseNodesPerDir = Get<Array<LO> >(fineLevel,
"lCoarseNodesPerDim");
145 Set<std::string>(coarseLevel,
"aggregationRegionType", regionType);
146 Set<int> (coarseLevel,
"numDimensions", numDimensions);
147 Set<Array<LO> > (coarseLevel,
"lNodesPerDim", lCoarseNodesPerDir);
149 if((pL.get<
bool>(
"interface aggregation") ==
true) && (regionType ==
"uncoupled")) {
150 Array<LO> coarseInterfacesDimensions = Get<Array<LO> >(fineLevel,
"coarseInterfacesDimensions");
151 Array<LO> nodeOnCoarseInterface = Get<Array<LO> >(fineLevel,
"nodeOnCoarseInterface");
152 Set<Array<LO> >(coarseLevel,
"interfacesDimensions", coarseInterfacesDimensions);
153 Set<Array<LO> >(coarseLevel,
"nodeOnInterface", nodeOnCoarseInterface);
156 }
else if(pL.get<
bool>(
"structured aggregation") ==
true) {
157 if(pL.get<
bool>(
"aggregation coupled") ==
true) {
158 gCoarseNodesPerDir = Get<Array<GO> >(fineLevel,
"gCoarseNodesPerDim");
159 Set<Array<GO> >(coarseLevel,
"gNodesPerDim", gCoarseNodesPerDir);
161 lCoarseNodesPerDir = Get<Array<LO> >(fineLevel,
"lCoarseNodesPerDim");
162 Set<Array<LO> >(coarseLevel,
"lNodesPerDim", lCoarseNodesPerDir);
163 numDimensions = Get<int>(fineLevel,
"numDimensions");
164 Set<int>(coarseLevel,
"numDimensions", numDimensions);
166 }
else if(pL.get<
bool>(
"Geometric") ==
true) {
167 coarseCoords = Get<RCP<xdMV> >(coarseLevel,
"coarseCoordinates");
168 gCoarseNodesPerDir = Get<Array<GO> >(coarseLevel,
"gCoarseNodesPerDim");
169 lCoarseNodesPerDir = Get<Array<LO> >(coarseLevel,
"lCoarseNodesPerDim");
170 Set<Array<GO> >(coarseLevel,
"gNodesPerDim", gCoarseNodesPerDir);
171 Set<Array<LO> >(coarseLevel,
"lNodesPerDim", lCoarseNodesPerDir);
173 Set<RCP<xdMV> >(coarseLevel,
"Coordinates", coarseCoords);
176 if (coarseLevel.
IsAvailable(
"Coordinates",
this)) {
177 GetOStream(
Runtime0) <<
"Reusing coordinates" << std::endl;
181 RCP<Aggregates> aggregates = Get< RCP<Aggregates> > (fineLevel,
"Aggregates");
182 fineCoords = Get< RCP<xdMV> >(fineLevel,
"Coordinates");
183 RCP<const Map> coarseMap = Get< RCP<const Map> > (fineLevel,
"CoarseMap");
189 ArrayView<const GO> elementAList = coarseMap->getNodeElementList();
192 if (rcp_dynamic_cast<const StridedMap>(coarseMap) != Teuchos::null)
193 blkSize = rcp_dynamic_cast<
const StridedMap>(coarseMap)->getFixedBlockSize();
195 GO indexBase = coarseMap->getIndexBase();
196 size_t numElements = elementAList.size() / blkSize;
197 Array<GO> elementList(numElements);
200 for (LO i = 0; i < Teuchos::as<LO>(numElements); i++)
201 elementList[i] = (elementAList[i*blkSize]-indexBase)/blkSize + indexBase;
203 RCP<const Map> uniqueMap = fineCoords->getMap();
204 RCP<const Map> coarseCoordMap = MapFactory ::Build(coarseMap->lib(), Teuchos::OrdinalTraits<Xpetra::global_size_t>::invalid(), elementList, indexBase, coarseMap->getComm());
205 coarseCoords = Xpetra::MultiVectorFactory<typename Teuchos::ScalarTraits<Scalar>::coordinateType,LO,GO,NO>::Build(coarseCoordMap, fineCoords->getNumVectors());
208 RCP<xdMV> ghostedCoords = fineCoords;
209 if (aggregates->AggregatesCrossProcessors()) {
210 RCP<const Map> nonUniqueMap = aggregates->GetMap();
211 RCP<const Import> importer = ImportFactory::Build(uniqueMap, nonUniqueMap);
213 ghostedCoords = Xpetra::MultiVectorFactory<typename Teuchos::ScalarTraits<Scalar>::coordinateType,LO,GO,NO>::Build(nonUniqueMap, fineCoords->getNumVectors());
214 ghostedCoords->doImport(*fineCoords, *importer, Xpetra::INSERT);
218 int myPID = uniqueMap->getComm()->getRank();
219 LO numAggs = aggregates->GetNumAggregates();
220 ArrayRCP<LO> aggSizes = aggregates->ComputeAggregateSizes();
221 const ArrayRCP<const LO> vertex2AggID = aggregates->GetVertex2AggId()->getData(0);
222 const ArrayRCP<const LO> procWinner = aggregates->GetProcWinner()->getData(0);
225 for (
size_t j = 0; j < fineCoords->getNumVectors(); j++) {
226 ArrayRCP<const typename Teuchos::ScalarTraits<Scalar>::coordinateType> fineCoordsData = ghostedCoords->getData(j);
227 ArrayRCP<typename Teuchos::ScalarTraits<Scalar>::coordinateType> coarseCoordsData = coarseCoords->getDataNonConst(j);
229 for (LO lnode = 0; lnode < vertex2AggID.size(); lnode++) {
230 if (procWinner[lnode] == myPID &&
231 lnode < vertex2AggID.size() &&
232 lnode < fineCoordsData.size() &&
233 vertex2AggID[lnode] < coarseCoordsData.size() &&
234 Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<Scalar>::coordinateType>::isnaninf(fineCoordsData[lnode]) ==
false) {
235 coarseCoordsData[vertex2AggID[lnode]] += fineCoordsData[lnode];
238 for (LO agg = 0; agg < numAggs; agg++) {
239 coarseCoordsData[agg] /= aggSizes[agg];
243 Set<RCP<xdMV> >(coarseLevel,
"Coordinates", coarseCoords);
246 int writeStart = pL.get<
int>(
"write start"), writeEnd = pL.get<
int>(
"write end");
247 if (writeStart == 0 && fineLevel.
GetLevelID() == 0 && writeStart <= writeEnd) {
248 std::ostringstream buf;
250 std::string fileName =
"coordinates_before_rebalance_level_" + buf.str() +
".m";
251 Xpetra::IO<typename Teuchos::ScalarTraits<Scalar>::coordinateType,LO,GO,NO>::Write(fileName,*fineCoords);
254 std::ostringstream buf;
256 std::string fileName =
"coordinates_before_rebalance_level_" + buf.str() +
".m";
257 Xpetra::IO<typename Teuchos::ScalarTraits<Scalar>::coordinateType,LO,GO,NO>::Write(fileName,*coarseCoords);
263 #endif // MUELU_COORDINATESTRANSFER_FACTORY_DEF_HPP RequestMode GetRequestMode() const
Timer to be used in factories. Similar to Monitor but with additional timers.
One-liner description of what is happening.
Namespace for MueLu classes and methods.
int GetLevelID() const
Return level number.
void DeclareInput(Level &finelevel, Level &coarseLevel) const
Specifies the data that this class needs, and the factories that generate that data.
Class that holds all level-specific information.
void Build(Level &fineLevel, Level &coarseLevel) const
Build an object with this factory.
bool IsAvailable(const std::string &ename, const FactoryBase *factory=NoFactory::get()) const
Test whether a need's value has been saved.
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.