10 #ifndef MUELU_CUTDROP_HPP
11 #define MUELU_CUTDROP_HPP
13 #include "Kokkos_Core.hpp"
14 #if KOKKOS_VERSION >= 40799
15 #include "KokkosKernels_ArithTraits.hpp"
17 #include "Kokkos_ArithTraits.hpp"
20 #include "MueLu_Utilities.hpp"
25 namespace MueLu::CutDrop {
37 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
52 #if KOKKOS_VERSION >= 40799
53 using ATS = KokkosKernels::ArithTraits<scalar_type>;
55 using ATS = Kokkos::ArithTraits<scalar_type>;
61 :
A(A_.getLocalMatrixDevice())
64 template <
class local_matrix_type2>
72 #if KOKKOS_VERSION >= 40799
73 using ATS = KokkosKernels::ArithTraits<scalar_type>;
75 using ATS = Kokkos::ArithTraits<scalar_type>;
79 const local_matrix_type2
A;
84 KOKKOS_INLINE_FUNCTION
87 ,
offset(A_.graph.row_map(rlid_))
90 KOKKOS_INLINE_FUNCTION
92 return ATS::magnitude(
A.values(
offset + x) *
A.values(
offset + x));
95 KOKKOS_INLINE_FUNCTION
118 KOKKOS_INLINE_FUNCTION
128 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node, Misc::StrengthMeasure measure>
137 using diag_view_type =
typename Kokkos::DualView<const scalar_type*, Kokkos::LayoutStride, typename Node::device_type, Kokkos::MemoryUnmanaged>::t_dev;
144 #if KOKKOS_VERSION >= 40799
145 using ATS = KokkosKernels::ArithTraits<scalar_type>;
147 using ATS = Kokkos::ArithTraits<scalar_type>;
156 :
A(A_.getLocalMatrixDevice())
160 auto lclDiag2d =
diagVec->getLocalViewDevice(Xpetra::Access::ReadOnly);
161 diag = Kokkos::subview(lclDiag2d, Kokkos::ALL(), 0);
164 auto lcl2d =
diagVec->getLocalViewDevice(Xpetra::Access::ReadOnly);
165 diag = Kokkos::subview(lcl2d, Kokkos::ALL(), 0);
169 template <
class local_matrix_type2,
class diag_view_type2>
177 #if KOKKOS_VERSION >= 40799
178 using ATS = KokkosKernels::ArithTraits<scalar_type>;
180 using ATS = Kokkos::ArithTraits<scalar_type>;
183 #if KOKKOS_VERSION >= 40799
184 using mATS = KokkosKernels::ArithTraits<magnitudeType>;
186 using mATS = Kokkos::ArithTraits<magnitudeType>;
189 const local_matrix_type2
A;
196 KOKKOS_INLINE_FUNCTION
201 ,
offset(A_.graph.row_map(rlid_))
204 KOKKOS_INLINE_FUNCTION
207 auto x_aij = ATS::magnitude(
A.values(
offset + x) *
A.values(
offset + x));
209 return (x_aij / x_aiiajj);
211 auto val =
A.values(
offset + x);
212 auto neg_aij = -ATS::real(val);
213 auto max_neg_aik = ATS::real(
diag(
rlid));
214 auto v = neg_aij / max_neg_aik;
215 if (ATS::real(v) <= mATS::zero()) {
216 return -ATS::magnitude(v * v);
218 return ATS::magnitude(v * v);
221 auto val =
A.values(
offset + x);
223 const bool is_nonpositive = ATS::real(val) <= mATS::zero();
224 magnitudeType aij2 = ATS::magnitude(val) * ATS::magnitude(val);
227 return (aij2 / x_aiiajj);
231 KOKKOS_INLINE_FUNCTION
254 KOKKOS_INLINE_FUNCTION
261 template <Misc::StrengthMeasure measure,
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
273 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node,
class DistanceFunctorType>
282 using diag_view_type =
typename Kokkos::DualView<const scalar_type*, Kokkos::LayoutStride, typename Node::device_type, Kokkos::MemoryUnmanaged>::t_dev;
289 #if KOKKOS_VERSION >= 40799
290 using ATS = KokkosKernels::ArithTraits<scalar_type>;
292 using ATS = Kokkos::ArithTraits<scalar_type>;
302 :
A(A_.getLocalMatrixDevice())
307 auto lclDiag2d =
diagVec->getLocalViewDevice(Xpetra::Access::ReadOnly);
308 diag = Kokkos::subview(lclDiag2d, Kokkos::ALL(), 0);
311 template <
class local_matrix_type2,
class DistanceFunctorType2,
class diag_view_type2>
319 #if KOKKOS_VERSION >= 40799
320 using ATS = KokkosKernels::ArithTraits<scalar_type>;
322 using ATS = Kokkos::ArithTraits<scalar_type>;
326 const local_matrix_type2
A;
336 KOKKOS_INLINE_FUNCTION
342 ,
offset(A_.graph.row_map(rlid_))
345 KOKKOS_INLINE_FUNCTION
347 auto clid =
A.graph.entries(
offset + x);
354 auto aij2 = ATS::magnitude(val) * ATS::magnitude(val);
358 KOKKOS_INLINE_FUNCTION
381 KOKKOS_INLINE_FUNCTION
391 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node,
class DistanceFunctorType, Misc::StrengthMeasure measure>
400 using diag_view_type =
typename Kokkos::DualView<const scalar_type*, Kokkos::LayoutStride, typename Node::device_type, Kokkos::MemoryUnmanaged>::t_dev;
407 #if KOKKOS_VERSION >= 40799
408 using ATS = KokkosKernels::ArithTraits<scalar_type>;
410 using ATS = Kokkos::ArithTraits<scalar_type>;
420 :
A(A_.getLocalMatrixDevice())
426 auto lclDiag2d =
diagVec->getLocalViewDevice(Xpetra::Access::ReadOnly);
427 diag = Kokkos::subview(lclDiag2d, Kokkos::ALL(), 0);
430 auto lclDiag2d =
diagVec->getLocalViewDevice(Xpetra::Access::ReadOnly);
431 diag = Kokkos::subview(lclDiag2d, Kokkos::ALL(), 0);
435 template <
class local_matrix_type2,
class DistanceFunctorType2,
class diag_view_type2>
443 #if KOKKOS_VERSION >= 40799
444 using ATS = KokkosKernels::ArithTraits<scalar_type>;
446 using ATS = Kokkos::ArithTraits<scalar_type>;
449 #if KOKKOS_VERSION >= 40799
450 using mATS = KokkosKernels::ArithTraits<magnitudeType>;
452 using mATS = Kokkos::ArithTraits<magnitudeType>;
455 const local_matrix_type2
A;
467 KOKKOS_INLINE_FUNCTION
473 ,
offset(A_.graph.row_map(rlid_))
476 KOKKOS_INLINE_FUNCTION
478 auto clid =
A.graph.entries(
offset + x);
487 auto aiiajj = ATS::magnitude(
diag(
rlid)) * ATS::magnitude(
diag(clid));
488 auto aij2 = ATS::magnitude(val) * ATS::magnitude(val);
489 return (aij2 / aiiajj);
491 auto neg_aij = -ATS::real(val);
492 auto max_neg_aik = ATS::real(
diag(
rlid));
493 auto v = ATS::magnitude(neg_aij / max_neg_aik);
494 if (ATS::real(neg_aij) >=
mzero)
499 auto aiiajj = ATS::magnitude(
diag(
rlid)) * ATS::magnitude(
diag(clid));
500 const bool is_nonpositive = ATS::real(val) <= mATS::zero();
501 magnitudeType aij2 = ATS::magnitude(val) * ATS::magnitude(val);
505 return aij2 / aiiajj;
509 KOKKOS_INLINE_FUNCTION
532 KOKKOS_INLINE_FUNCTION
539 template <Misc::StrengthMeasure measure,
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node,
class DistanceFunctorType>
541 DistanceFunctorType& dist2_,
556 template <
class comparison_type>
565 #if KOKKOS_VERSION >= 40799
566 using ATS = KokkosKernels::ArithTraits<scalar_type>;
568 using ATS = Kokkos::ArithTraits<scalar_type>;
577 Kokkos::View<local_ordinal_type*, memory_space>
index;
585 index = Kokkos::View<local_ordinal_type*, memory_space>(
"indices",
A.nnz());
588 KOKKOS_FORCEINLINE_FUNCTION
590 auto row =
A.rowConst(rlid);
591 size_t nnz = row.length;
593 auto drop_view = Kokkos::subview(
results, Kokkos::make_pair(
A.graph.row_map(rlid),
A.graph.row_map(rlid + 1)));
594 auto row_permutation = Kokkos::subview(
index, Kokkos::make_pair(
A.graph.row_map(rlid),
A.graph.row_map(rlid + 1)));
596 auto comparator =
comparison.getComparator(rlid);
598 #ifdef MUELU_COALESCE_DROP_DEBUG
600 Kokkos::printf(
"SoC: ");
602 Kokkos::printf(
"%5f ", comparator.get_value(k));
604 Kokkos::printf(
"\n");
608 for (
size_t i = 0; i < nnz; ++i) {
609 row_permutation(i) = i;
613 size_t keepStart = 0;
614 size_t dropStart = nnz;
616 for (
size_t i = 1; i < nnz; ++i) {
617 auto const& x = row_permutation(i - 1);
618 auto const& y = row_permutation(i);
625 if (
eps *
eps * x_aij > y_aij) {
633 for (
size_t i = keepStart; i < nnz; ++i) {
634 drop_view(row_permutation(i)) = Kokkos::max(dropStart <= i ?
DROP :
KEEP, drop_view(row_permutation(i)));
Kokkos::View< const bool *, memory_space > boundary_nodes_view
KOKKOS_INLINE_FUNCTION bool operator()(size_t x, size_t y) const
const diag_view_type2 diag
typename local_matrix_type::memory_space memory_space
auto make_comparison_functor(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A_, typename ScaledComparison< Scalar, LocalOrdinal, GlobalOrdinal, Node, measure >::results_view &results_)
static Teuchos::RCP< Vector > GetMatrixMaxMinusOffDiagonal(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A)
Return vector containing: max_{i=k}(-a_ik), for each for i in the matrix.
const local_ordinal_type offset
typename local_matrix_type::ordinal_type local_ordinal_type
UnscaledComparison(matrix_type &A_, results_view &results_)
const results_view results
Comparator< local_matrix_type, diag_view_type > comparator_type
ScaledDistanceLaplacianComparison(matrix_type &A_, DistanceFunctorType &dist2_, results_view &results_)
const local_matrix_type2 A
KOKKOS_INLINE_FUNCTION magnitudeType get_value(size_t x) const
typename local_matrix_type2::ordinal_type local_ordinal_type
KOKKOS_INLINE_FUNCTION comparator_type getComparator(local_ordinal_type rlid) const
const local_matrix_type2 A
typename ATS::magnitudeType magnitudeType
KOKKOS_INLINE_FUNCTION comparator_type getComparator(local_ordinal_type rlid) const
Kokkos::ArithTraits< magnitudeType > mATS
Kokkos::View< DecisionType *, memory_space > results_view
typename local_matrix_type2::memory_space memory_space
Kokkos::View< DecisionType *, memory_space > results_view
typename matrix_type::local_matrix_type local_matrix_type
Teuchos::RCP< diag_vec_type > diagVec
KOKKOS_INLINE_FUNCTION Comparator(const local_matrix_type2 &A_, const diag_view_type2 &diag_, const DistanceFunctorType2 *dist2_, local_ordinal_type rlid_, const results_view &results_)
Kokkos::View< local_ordinal_type *, memory_space > index
Kokkos::View< DecisionType *, memory_space > results_view
typename local_matrix_type::ordinal_type local_ordinal_type
typename local_matrix_type::ordinal_type local_ordinal_type
Teuchos::RCP< diag_vec_type > diagVec
Teuchos::RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getMaxMinusOffDiagonal(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, DistanceFunctorType &distFunctor)
const diag_view_type2 diag
typename Kokkos::DualView< const scalar_type *, Kokkos::LayoutStride, typename Node::device_type, Kokkos::MemoryUnmanaged >::t_dev diag_view_type
typename local_matrix_type2::value_type scalar_type
Kokkos::View< DecisionType *, memory_space > results_view
Orders entries of row by .
const local_matrix_type2 A
Kokkos::ArithTraits< scalar_type > ATS
typename local_matrix_type::memory_space memory_space
ScaledComparison(matrix_type &A_, results_view &results_)
KOKKOS_INLINE_FUNCTION comparator_type getComparator(local_ordinal_type rlid) const
comparison_type comparison
const local_ordinal_type rlid
Orders entries of row by where is the distance Laplacian.
Comparator< local_matrix_type > comparator_type
typename local_matrix_type::ordinal_type local_ordinal_type
typename ATS::magnitudeType magnitudeType
typename matrix_type::local_matrix_type local_matrix_type
typename comparison_type::local_matrix_type local_matrix_type
const local_ordinal_type offset
const local_ordinal_type rlid
const local_ordinal_type rlid
Kokkos::View< DecisionType *, memory_space > results_view
Comparator< local_matrix_type, DistanceFunctorType, diag_view_type > comparator_type
typename local_matrix_type::value_type scalar_type
typename local_matrix_type::value_type scalar_type
typename local_matrix_type::ordinal_type local_ordinal_type
KOKKOS_INLINE_FUNCTION void serialHeapSort(view_type &v, comparator_type comparator)
typename local_matrix_type2::ordinal_type local_ordinal_type
KOKKOS_INLINE_FUNCTION bool operator()(size_t x, size_t y) const
Kokkos::ArithTraits< scalar_type > ATS
const local_matrix_type2 A
Kokkos::ArithTraits< magnitudeType > mATS
Kokkos::ArithTraits< scalar_type > ATS
typename local_matrix_type2::ordinal_type local_ordinal_type
UnscaledDistanceLaplacianComparison(matrix_type &A_, DistanceFunctorType &dist2_, results_view &results_)
typename Kokkos::DualView< const scalar_type *, Kokkos::LayoutStride, typename Node::device_type, Kokkos::MemoryUnmanaged >::t_dev diag_view_type
Kokkos::View< DecisionType *, memory_space > results_view
Kokkos::ArithTraits< scalar_type > ATS
typename Kokkos::DualView< const scalar_type *, Kokkos::LayoutStride, typename Node::device_type, Kokkos::MemoryUnmanaged >::t_dev diag_view_type
Kokkos::View< DecisionType *, memory_space > results_view
Kokkos::ArithTraits< scalar_type > ATS
typename ATS::magnitudeType magnitudeType
Order each row by a criterion, compare the ratio of values and drop all entries once the ratio is bel...
KOKKOS_INLINE_FUNCTION magnitudeType get_value(size_t x) const
typename ATS::magnitudeType magnitudeType
typename local_matrix_type::memory_space memory_space
const results_view results
typename ATS::magnitudeType magnitudeType
typename local_matrix_type::value_type scalar_type
static RCP< Vector > GetMatrixOverlappedDiagonal(const Matrix &A)
Extract Overlapped Matrix Diagonal.
Kokkos::ArithTraits< scalar_type > ATS
auto make_dlap_comparison_functor(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A_, DistanceFunctorType &dist2_, typename ScaledDistanceLaplacianComparison< Scalar, LocalOrdinal, GlobalOrdinal, Node, DistanceFunctorType, measure >::results_view &results_)
typename matrix_type::local_matrix_type local_matrix_type
const results_view results
Teuchos::RCP< diag_vec_type > diagVec
const DistanceFunctorType2 * dist2
KOKKOS_FORCEINLINE_FUNCTION void operator()(const local_ordinal_type &rlid) const
typename local_matrix_type::value_type scalar_type
typename local_matrix_type2::value_type scalar_type
KOKKOS_INLINE_FUNCTION comparator_type getComparator(local_ordinal_type rlid) const
typename ATS::magnitudeType magnitudeType
KOKKOS_INLINE_FUNCTION bool operator()(size_t x, size_t y) const
typename local_matrix_type2::ordinal_type local_ordinal_type
typename local_matrix_type2::memory_space memory_space
DistanceFunctorType dist2
Kokkos::ArithTraits< scalar_type > ATS
typename ATS::magnitudeType magnitudeType
Kokkos::ArithTraits< scalar_type > ATS
KOKKOS_INLINE_FUNCTION bool operator()(size_t x, size_t y) const
typename local_matrix_type2::memory_space memory_space
KOKKOS_INLINE_FUNCTION magnitudeType get_value(size_t x) const
const diag_view_type2 diag
const results_view results
const local_ordinal_type offset
Kokkos::ArithTraits< scalar_type > ATS
const DistanceFunctorType2 * dist2
typename local_matrix_type2::value_type scalar_type
Teuchos::RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getDiagonal(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, DistanceFunctorType &distFunctor)
KOKKOS_INLINE_FUNCTION Comparator(const local_matrix_type2 &A_, const diag_view_type2 &diag_, const local_ordinal_type rlid_, const results_view &results_)
KOKKOS_INLINE_FUNCTION Comparator(const local_matrix_type2 &A_, local_ordinal_type rlid_, const results_view &results_)
typename ATS::magnitudeType magnitudeType
Comparator< local_matrix_type, DistanceFunctorType, diag_view_type > comparator_type
KOKKOS_INLINE_FUNCTION magnitudeType get_value(size_t x) const
CutDropFunctor(comparison_type &comparison_, magnitudeType threshold)
typename ATS::magnitudeType magnitudeType
typename matrix_type::local_matrix_type local_matrix_type
Kokkos::View< DecisionType *, memory_space > results_view
typename local_matrix_type::memory_space memory_space
Orders entries of row by .
const local_ordinal_type offset
const magnitudeType mzero
typename local_matrix_type::memory_space memory_space
Kokkos::View< DecisionType *, memory_space > results_view
KOKKOS_INLINE_FUNCTION Comparator(const local_matrix_type2 &A_, const diag_view_type2 &diag_, const DistanceFunctorType2 *dist2_, local_ordinal_type rlid_, const results_view &results_)
typename local_matrix_type::value_type scalar_type
typename local_matrix_type2::value_type scalar_type
typename local_matrix_type2::memory_space memory_space
DistanceFunctorType dist2