46#ifndef MUELU_COORDINATESTRANSFER_FACTORY_KOKKOS_DEF_HPP
47#define MUELU_COORDINATESTRANSFER_FACTORY_KOKKOS_DEF_HPP
51#include <Xpetra_ImportFactory.hpp>
52#include <Xpetra_IO.hpp>
53#include <Xpetra_MultiVectorFactory.hpp>
54#include <Xpetra_MapFactory.hpp>
56#include "MueLu_Aggregates_kokkos.hpp"
57#include "MueLu_CoarseMapFactory_kokkos.hpp"
60#include "MueLu_Utilities_kokkos.hpp"
64 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class DeviceType>
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<
int > (
"write start", -1,
"First level at which coordinates should be written to file");
72 validParamList->set<
int > (
"write end", -1,
"Last level at which coordinates should be written to file");
73 validParamList->set<
bool > (
"structured aggregation",
false,
"Flag specifying that the geometric data is transferred for StructuredAggregationFactory");
74 validParamList->set<RCP<const FactoryBase> > (
"lCoarseNodesPerDim", Teuchos::null,
"Factory providing the local number of nodes per spatial dimensions of the mesh");
75 validParamList->set<RCP<const FactoryBase> > (
"numDimensions", Teuchos::null,
"Factory providing the number of spatial dimensions of the mesh");
77 return validParamList;
80 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class DeviceType>
82 static bool isAvailableCoords =
false;
84 const ParameterList& pL = GetParameterList();
85 if(pL.get<
bool>(
"structured aggregation") ==
true) {
86 Input(fineLevel,
"lCoarseNodesPerDim");
87 Input(fineLevel,
"numDimensions");
90 isAvailableCoords = coarseLevel.
IsAvailable(
"Coordinates",
this);
92 if (isAvailableCoords ==
false) {
93 Input(fineLevel,
"Coordinates");
94 Input(fineLevel,
"Aggregates");
95 Input(fineLevel,
"CoarseMap");
100 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class DeviceType>
104 typedef Xpetra::MultiVector<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LO,GO,NO> doubleMultiVector;
106 GetOStream(
Runtime0) <<
"Transferring coordinates" << std::endl;
108 if (coarseLevel.
IsAvailable(
"Coordinates",
this)) {
109 GetOStream(
Runtime0) <<
"Reusing coordinates" << std::endl;
114 const ParameterList& pL = GetParameterList();
115 if(pL.get<
bool>(
"structured aggregation") ==
true) {
116 Array<LO> lCoarseNodesPerDir = Get<Array<LO> >(fineLevel,
"lCoarseNodesPerDim");
117 Set<Array<LO> >(coarseLevel,
"lNodesPerDim", lCoarseNodesPerDir);
118 int numDimensions = Get<int>(fineLevel,
"numDimensions");
119 Set<int>(coarseLevel,
"numDimensions", numDimensions);
124 auto aggregates = Get<RCP<Aggregates_kokkos> >(fineLevel,
"Aggregates");
125 auto fineCoords = Get<RCP<doubleMultiVector> >(fineLevel,
"Coordinates");
126 auto coarseMap = Get<RCP<const Map> > (fineLevel,
"CoarseMap");
132 ArrayView<const GO> elementAList = coarseMap->getLocalElementList();
133 GO indexBase = coarseMap->getIndexBase();
136 if (rcp_dynamic_cast<const StridedMap>(coarseMap) != Teuchos::null)
137 blkSize = rcp_dynamic_cast<const StridedMap>(coarseMap)->getFixedBlockSize();
139 Array<GO> elementList;
140 ArrayView<const GO> elementListView;
144 elementListView = elementAList;
147 auto numElements = elementAList.size() / blkSize;
149 elementList.resize(numElements);
152 for (LO i = 0; i < Teuchos::as<LO>(numElements); i++)
153 elementList[i] = (elementAList[i*blkSize]-indexBase)/blkSize + indexBase;
155 elementListView = elementList;
158 auto uniqueMap = fineCoords->getMap();
159 auto coarseCoordMap = MapFactory::Build(coarseMap->lib(), Teuchos::OrdinalTraits<Xpetra::global_size_t>::invalid(),
160 elementListView, indexBase, coarseMap->getComm());
161 RCP<doubleMultiVector> coarseCoords = Xpetra::MultiVectorFactory<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LO,GO,NO>::Build(coarseCoordMap, fineCoords->getNumVectors());
164 RCP<doubleMultiVector> ghostedCoords = fineCoords;
165 if (aggregates->AggregatesCrossProcessors()) {
166 auto nonUniqueMap = aggregates->GetMap();
167 auto importer = ImportFactory::Build(uniqueMap, nonUniqueMap);
169 ghostedCoords = Xpetra::MultiVectorFactory<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LO,GO,NO>::Build(nonUniqueMap, fineCoords->getNumVectors());
170 ghostedCoords->doImport(*fineCoords, *importer, Xpetra::INSERT);
175 auto aggGraph = aggregates->GetGraph();
176 auto numAggs = aggGraph.numRows();
178 auto fineCoordsView = fineCoords ->getDeviceLocalView(Xpetra::Access::ReadOnly);
179 auto coarseCoordsView = coarseCoords->getDeviceLocalView(Xpetra::Access::OverwriteAll);
185 const auto dim = fineCoords->getNumVectors();
187 typename AppendTrait<
decltype(fineCoordsView), Kokkos::RandomAccess>::type fineCoordsRandomView = fineCoordsView;
188 for (
size_t j = 0; j < dim; j++) {
189 Kokkos::parallel_for(
"MueLu:CoordinatesTransferF:Build:coord", Kokkos::RangePolicy<local_ordinal_type, execution_space>(0, numAggs),
190 KOKKOS_LAMBDA(
const LO i) {
194 auto aggregate = aggGraph.rowConst(i);
196 typename Teuchos::ScalarTraits<Scalar>::magnitudeType sum = 0.0;
197 for (
size_t colID = 0; colID < static_cast<size_t>(aggregate.length); colID++)
198 sum += fineCoordsRandomView(aggregate(colID),j);
200 coarseCoordsView(i,j) = sum / aggregate.length;
205 Set<RCP<doubleMultiVector> >(coarseLevel,
"Coordinates", coarseCoords);
207 int writeStart = pL.get<
int>(
"write start"), writeEnd = pL.get<
int>(
"write end");
208 if (writeStart == 0 && fineLevel.
GetLevelID() == 0 && writeStart <= writeEnd) {
209 std::string fileName =
"coordinates_before_rebalance_level_" +
toString(fineLevel.
GetLevelID()) +
".m";
210 Xpetra::IO<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LO,GO,NO>::Write(fileName, *fineCoords);
213 std::string fileName =
"coordinates_before_rebalance_level_" +
toString(coarseLevel.
GetLevelID()) +
".m";
214 Xpetra::IO<typename Teuchos::ScalarTraits<Scalar>::magnitudeType,LO,GO,NO>::Write(fileName,*coarseCoords);
Timer to be used in factories. Similar to Monitor but with additional timers.
Class that holds all level-specific information.
bool IsAvailable(const std::string &ename, const FactoryBase *factory=NoFactory::get()) const
Test whether a need's value has been saved.
int GetLevelID() const
Return level number.
RequestMode GetRequestMode() const
Timer to be used in factories. Similar to SubMonitor but adds a timer level by level.
Namespace for MueLu classes and methods.
@ Runtime0
One-liner description of what is happening.
std::string toString(const T &what)
Little helper function to convert non-string types to strings.