46 #ifndef MUELU_MATRIXFREETENTATIVEP_DEF_HPP
47 #define MUELU_MATRIXFREETENTATIVEP_DEF_HPP
51 #include "MueLu_Aggregates.hpp"
53 #include <Tpetra_KokkosCompat_ClassicNodeAPI_Wrapper.hpp>
57 #include "Xpetra_Map.hpp"
65 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
71 using impl_scalar_type =
typename Kokkos::ArithTraits<Scalar>::val_type;
72 impl_scalar_type implAlpha = alpha;
78 typename Aggregates::aggregates_sizes_type::const_type aggSizes = aggregates_->ComputeAggregateSizes();
80 auto kokkos_view_X = X.getDeviceLocalView(Xpetra::Access::ReadOnly);
81 auto kokkos_view_Y = Y.getDeviceLocalView(Xpetra::Access::ReadWrite);
82 LO numCols = kokkos_view_X.extent(1);
85 auto vertex2AggId = aggregates_->GetVertex2AggId();
86 auto vertex2AggIdView = vertex2AggId->getDeviceLocalView(Xpetra::Access::ReadOnly);
87 LO numNodes = kokkos_view_X.extent(0);
92 "MueLu:MatrixFreeTentativeR_kokkos:apply",
md_range_type({0, 0}, {numCols, numNodes}),
93 KOKKOS_LAMBDA(
const int colIdx,
const int NodeIdx) {
94 LO aggIdx = vertex2AggIdView(NodeIdx, 0);
96 Kokkos::atomic_add(&kokkos_view_Y(aggIdx, colIdx), implAlpha * kokkos_view_X(NodeIdx, colIdx) / Kokkos::sqrt(aggSizes(aggIdx)));
100 const auto vertex2Agg = aggregates_->GetVertex2AggId();
101 auto vertex2AggView = vertex2Agg->getDeviceLocalView(Xpetra::Access::ReadOnly);
102 LO numNodes = kokkos_view_Y.extent(0);
106 Kokkos::parallel_for(
107 "MueLu:MatrixFreeTentativeP:apply",
md_range_type({0, 0}, {numCols, numNodes}),
108 KOKKOS_LAMBDA(
const int colIdx,
const int fineIdx) {
109 LO aggIdx = vertex2AggView(fineIdx, 0);
110 kokkos_view_Y(fineIdx, colIdx) += implAlpha * kokkos_view_X(aggIdx, colIdx) / Kokkos::sqrt(aggSizes(aggIdx));
116 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
123 #define MUELU_MATRIXFREETENTATIVEP_SHORT
124 #endif // MUELU_MATRIXFREETENTATIVEP_DEF_HPP
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void residual(const MultiVector &X, const MultiVector &B, MultiVector &R) const override
MueLu::DefaultScalar Scalar
virtual void scale(const Scalar &alpha)=0
Kokkos::MDRangePolicy< local_ordinal_type, execution_space, Kokkos::Rank< 2 > > md_range_type
Exception throws to report errors in the internal logical of the program.
void apply(const MultiVector &X, MultiVector &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=Teuchos::ScalarTraits< Scalar >::one(), Scalar beta=Teuchos::ScalarTraits< Scalar >::zero()) const override