47 #ifndef MUELU_AGGREGATES_KOKKOS_DEF_HPP 48 #define MUELU_AGGREGATES_KOKKOS_DEF_HPP 50 #include <Xpetra_Map.hpp> 51 #include <Xpetra_Vector.hpp> 52 #include <Xpetra_VectorFactory.hpp> 54 #include "MueLu_LWGraph_kokkos.hpp" 61 template <
class LocalOrdinal,
class GlobalOrdinal,
class DeviceType>
62 Aggregates_kokkos<LocalOrdinal, GlobalOrdinal, Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> >::Aggregates_kokkos(LWGraph_kokkos graph) {
65 vertex2AggId_ = LOVectorFactory::Build(graph.GetImportMap());
68 procWinner_ = LOVectorFactory::Build(graph.GetImportMap());
71 isRoot_ = Teuchos::ArrayRCP<bool>(graph.GetImportMap()->getNodeNumElements(),
false);
74 aggregatesIncludeGhosts_ =
true;
78 template <
class LocalOrdinal,
class GlobalOrdinal,
class DeviceType>
79 Aggregates_kokkos<LocalOrdinal, GlobalOrdinal, Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> >::Aggregates_kokkos(
const RCP<const Map>& map) {
82 vertex2AggId_ = LOVectorFactory::Build(map);
85 procWinner_ = LOVectorFactory::Build(map);
88 isRoot_ = Teuchos::ArrayRCP<bool>(map->getNodeNumElements(),
false);
91 aggregatesIncludeGhosts_ =
true;
95 template <
class LocalOrdinal,
class GlobalOrdinal,
class DeviceType>
96 typename Aggregates_kokkos<LocalOrdinal, GlobalOrdinal, Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> >::aggregates_sizes_type::const_type
97 Aggregates_kokkos<LocalOrdinal, GlobalOrdinal, Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> >::ComputeAggregateSizes(
bool forceRecompute,
bool cacheSizes)
const {
98 if (aggregateSizes_.size() && !forceRecompute) {
99 return aggregateSizes_;
104 aggregateSizes_ = aggregates_sizes_type(
"aggregates", 0);
106 aggregates_sizes_type aggregateSizes(
"aggregates", nAggregates_);
108 int myPID = GetMap()->getComm()->getRank();
110 auto vertex2AggId = vertex2AggId_->template getLocalView<DeviceType>();
111 auto procWinner = procWinner_ ->template getLocalView<DeviceType>();
113 typename AppendTrait<decltype(aggregateSizes_), Kokkos::Atomic>::type aggregateSizesAtomic = aggregateSizes;
114 Kokkos::parallel_for(
"Aggregates:ComputeAggregateSizes:for", procWinner.size(), KOKKOS_LAMBDA(
const LO k) {
115 if (procWinner(k, 0) == myPID)
116 aggregateSizesAtomic(vertex2AggId(k, 0))++;
120 aggregateSizes_ = aggregateSizes;
122 return aggregateSizes;
127 template <
class LocalOrdinal,
class GlobalOrdinal,
class DeviceType>
128 typename Aggregates_kokkos<LocalOrdinal, GlobalOrdinal, Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> >::local_graph_type
129 Aggregates_kokkos<LocalOrdinal, GlobalOrdinal, Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> >::GetGraph()
const {
130 typedef typename local_graph_type::row_map_type row_map_type;
131 typedef typename local_graph_type::entries_type entries_type;
133 int myPID = GetMap()->getComm()->getRank();
135 ArrayRCP<LO> vertex2AggId = vertex2AggId_->getDataNonConst(0);
136 ArrayRCP<LO> procWinner = procWinner_->getDataNonConst(0);
138 typename aggregates_sizes_type::const_type sizes = ComputeAggregateSizes();
140 int numAggregates = nAggregates_;
142 typename row_map_type::non_const_type rows(
"row_map", numAggregates+1);
143 for (LO i = 0; i < nAggregates_; i++)
144 rows(i+1) = rows(i) + sizes(i);
146 aggregates_sizes_type offsets(
"offsets", numAggregates);
147 for (LO i = 0; i < numAggregates; i++)
148 offsets(i) = rows(i);
150 typename entries_type::non_const_type cols(
"entries", rows(nAggregates_));
151 for (LO i = 0; i < procWinner.size(); i++)
152 if (procWinner[i] == myPID)
153 cols(offsets(vertex2AggId[i])++) = i;
155 return local_graph_type(cols, rows);
158 template <
class LocalOrdinal,
class GlobalOrdinal,
class DeviceType>
159 std::string Aggregates_kokkos<LocalOrdinal, GlobalOrdinal, Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> >::description()
const {
163 template <
class LocalOrdinal,
class GlobalOrdinal,
class DeviceType>
164 void Aggregates_kokkos<LocalOrdinal, GlobalOrdinal, Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> >::print(Teuchos::FancyOStream& out,
const Teuchos::EVerbosityLevel verbLevel)
const {
168 out0 <<
"Global number of aggregates: " << GetNumGlobalAggregates() << std::endl;
171 template <
class LocalOrdinal,
class GlobalOrdinal,
class DeviceType>
172 GlobalOrdinal Aggregates_kokkos<LocalOrdinal, GlobalOrdinal, Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType> >::GetNumGlobalAggregates()
const {
173 LO nAggregates = GetNumAggregates();
174 GO nGlobalAggregates;
175 MueLu_sumAll(vertex2AggId_->getMap()->getComm(), (GO)nAggregates, nGlobalAggregates);
176 return nGlobalAggregates;
179 template <
class LocalOrdinal,
class GlobalOrdinal,
class DeviceType>
180 const RCP<const Xpetra::Map<LocalOrdinal,GlobalOrdinal, Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType>> >
181 Aggregates_kokkos<LocalOrdinal, GlobalOrdinal, Kokkos::Compat::KokkosDeviceWrapperNode<DeviceType>>::GetMap()
const {
182 return vertex2AggId_->getMap();
187 #endif // MUELU_AGGREGATES_KOKKOS_DEF_HPP
#define MueLu_sumAll(rcpComm, in, out)
std::string toString(const T &what)
Little helper function to convert non-string types to strings.
Namespace for MueLu classes and methods.
Print statistics that do not involve significant additional computation.
#define MUELU_UNAGGREGATED
#define MUELU_DESCRIBE
Helper macro for implementing Describable::describe() for BaseClass objects.
void parallel_for(const ExecPolicy &policy, const FunctorType &functor, const std::string &str="", typename Impl::enable_if< !Impl::is_integral< ExecPolicy >::value >::type *=0)
virtual std::string description() const
Return a simple one-line description of this object.